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Abstract. The electronic and magnetic properties of many strongly-correlated 
systems are controlled by a limited number of states, located near the Fermi level 
and well isolated from the rest of the spectrum. This opens a formal way for 
combining the methods of first-principles electronic structure calculations, based on the 
C^ ' density-functional theory (DFT), with many-body models, formulated in the restricted 

S , Hilbert space of states close to the Fermi level. The core of this project is the SO- 

IL ' called "realistic modeling" or the construction of the model many-body Hamiltonians 

(-H . entirely from the first principles. Such a construction should be able to go beyond 

O I the conventional local-density approximation (LDA), which typically supplements the 

I ^i ' density- functional theory, and incorporate the physics of Coulomb correlations. It 

should also provide a transparent physical picture for the low-energy properties of 
strongly correlated materials. In this review article, we will outline the basic ideas of 
such a realistic modeling, which consists of the following steps: (i) The construction 
[~^ ] of the complete Wannier basis set for the low-energy LDA band; (ii) The construction 

Cn ' of the one-electron part of the model Hamiltonian in this Wannier basis; (iii) The 

^ . . calculation of the screened Coulomb interactions for the low-energy bands by means 

r^ I of the constrained DFT. The most difficult part of this project is the evaluation of 

>— s ' the screening caused by outer bands, which may have the same (e.g., the transition- 

metal 3d) character as the low-energy bands. The latter part can be efficiently done by 
combining the constrained DFT with the random-phase approximation for the screened 
Coulomb interaction. The entire procedure will be illustrated on the series of examples, 
including the distorted transition-metal perovskite oxides, the compounds with the 
C^ ' inversion symmetry breaking caused by the defects, and the alkali hyperoxide KO2, 

which can be regarded as an analog of strongly-correlated systems where the localized 
electrons reside on the molecular orbitals of the Oj^ dimer. In order to illustrate 
abilities of the realistic modeling, we will also consider solutions of the obtained low- 
energy models for a number of systems, and argue that it can be used as a powerful tool 
for the exploration and understanding of properties of strongly correlated materials. 



PACS numbers: 71.15.-m, 71.28. -Fd, 71.10.-w, 75.25. -Fz 
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1. Introduction 

Many successes of modern condensed-matter physics and cliemistry are related with 
the development of the density- functional theory (DFT), which is designed for the 
exploration of the ground state properties of various substances and based on the 
minimization of the total energy functional E[p\ with respect to the electron density 
p [H [21 [3] . For practical applications, DFT resorts to iterative solution of one-electron 
Kohn-Sham equations 

-|^V^ + V^i;,= 5,^„ (1) 

together with the equation for the electron density 

P = J2fM^ (2) 

i 

defined in terms of eigenfunctions (ipi), eigenvalues (£j), and occupation numbers (/j) of 
Kohn-Sham quasiparticles. The potential V can be divided into the Coulomb (Vh), 
exchange- correlation (Vxc), and the external parts (14xt), which are the functional 
derivatives of corresponding contributions to the total energy with respect to the electron 
density. Formally speaking, this procedure is fully ab initio and free of any adjustable 
parameters. 

However, the form of the exchange- correlation potential is generally unknown. For 
practical purposes, it is typically treated in the local-density approximation (LDA), 
which employs an analytical expression borrowed from the theory of homogeneous 
electron gas in which the density of the electron gas is replaced by the local density 
of the real system. LDA is far from being perfect and there are many examples of the 
so-called strongly correlated materials where the conventional LDA fails in describing 
the excited- as well as the ground-state properties ^. 

In the strongly correlated materials, the state of each electron strongly depends 
on the state of other electrons of the system, which are coupled (or correlate with 
each other) via the Coulomb interaction. Thus, this is the real many-body problem, 
and the situation is very different from the behavior of the homogeneous electron gas. 
The canonical example of strongly correlated materials is the transition-metal oxides 
[1] . A typical example of the electronic structure of the transition-metal oxides in the 
local-density approximation is shown in Figure [1] for the series of distorted perovskite 
compounds!^ We would like to emphasize two points. 

(i) The common feature of many strongly correlated systems is the existence of a 
limited group of states, located near the Fermi level and well isolated from the 
rest of the spectrum. In the following, these states will be called the as "low-energy 
states" or the "low-energy bands" or, simply, the L-bands. In the case of perovskite 
oxides depicted in Figure [1] these are the narrow transition-metal t2g bands. From 
this point of view, the theoretical description of the strongly correlated systems 



The properties of the distorted perovskite oxides wih be discussed in details in Section 16.3 
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Figure 1. Total and partial densities of states of YTiOa (a), LaTiOs (b), YVO3 
(c, orthorhombic phase), and LaVOs (d) in the local-density approximation |S]. The 
shaded area shows the contributions of transition-metal 3d-states. Other symbols show 
the positions of the main bands. The Fermi level is at zero energy. 



is feasible, and this is certainly a good sign. For example, if we are interested in 
electronic or magnetic properties, that are mainly controlled by the states close 
to the Fermi level, we can mainly concentrate on the behavior of this group of 
states and disengage ourself from other details of the electronic structure. Since 
the number of such states is limited, the problem can be solved, at least numerically. 

(ii) However, the bad points is that in order to solve this problem we should inevitably 
go beyond the local-density approximation, which greatly oversimplifies the physics 
of Coulomb correlations. For example, the systems depicted in Figure [H are metals 
within LDA, while in practice all of them are Mott insulators |1]. 

The insulating behavior is frequently associated with the excited state properties, 
which are not supposed to be reproduced by the Kohn-Sham equations designed for 
the ground state. However, the problem is much more serious. Suppose that we 
are interested in the behavior of interatomic magnetic interactions, which are the 
ground state properties. For the narrow-band compounds, the main contributions to 
the magnetic interactions can be typically identified with the double exchange and 
superexchange mechanisms [6l[7j. The double exchange operates in the metallic systems. 
It is always ferromagnetic and proportional to the kinetic-energy gain, t, caused by free 
electron hoppings between different sites of the lattice [8] . The superexchange dominates 
in insulating compounds, where the double exchange vanishes [9]. The superexchange 
can be either antiferromagnetic or ferromagnetic, depending on the number of electrons 
and the orbital state of the system ^U\ . It is proportional to t^ and inversely proportional 
to the parameter of on-site Coulomb interaction U [llj. Now, suppose that because 
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of the limitations of LDA, our system is metallic rather than insulating. Then, the 
ferromagnetic double exchange will clearly dominate and we may have a totally wrong 
answer for the interatomic magnetic interactions as well as for the magnetic ground 
state. Therefore, if the system is insulating, an energy gap should be present already 
in the spectrum of Kohn-Sham eigenvalues. If it is not, the approximations underlying 
these Kohn-Sham equations must be corrected. 

Because of these difficulties, for the long time the strongly correlated compounds 
have been almost exclusive prerogative of the model condensed-matter physics, which 
has accumulated a great knowledge about treating and solution of this many-body 
problem but only for model systems, described in terms of a limited number of model 
parameters. The model analysis has indeed provided a useful and insightful information 
about many properties of the strongly correlated systems. However, as the complexity 
of such systems rose, it inevitably posed a question about the correct choice of the 
model parameters, and in a number of cases - of the model itself. For example, the 
typical questions are how to incorporate the information about the chemical signature of 
elements composing the lattice into the model or how to treat lattice distortions? These 
questions gave rise to the new direction, which can be called as the "realistic modeling" . 
The basic idea of the realistic modeling is to construct the model Hamiltonian entirely 
from the first principles, and make it free of any adjustable parameters, and then to 
solve this model Hamiltonian by using modern many-body methods. This was the main 
motivation in the area of model condensed-matter physics. 

Then, why do we need the realistic modeling in the area of computational physics 
dealing with the first-principles electronic structure calculations? 

(i) Presumably, it is needless to say that the solution of the many-body problem for 
many realistic systems is a tremendous computational task. Therefore, one would 
always like to formulate it rigorously only in a restricted Hilbert space, which 
would pick up only those states that are primarily responsible for the physics of 
either of the considered phenomenon, and include the effect of other states into the 
renormalization of the model parameters in the primary Hilbert space. 

(ii) It is true that the methods of first-principles electronic structure calculations are 
currently on the rise. The standard computations within the density-functional 
theory become a matter of routine and frequently regarded as a useful tool for the 
numerical modeling or the numerical experiment in the materials science. However, 
we will always need some additional tools for the analysis and interpretation of 
the obtained data. Our final goal is not to obtain or reproduce the experimental 
trend for some complex compounds. Our goal is provide some transparent physical 
interpretation for this trend on the microscopic level and come up with some useful 
suggestions about how this trend can be further improved. In this sense, the realistic 
modeling will continue to play an important role as the tool for the analysis and 
interpretation of data obtained in the standard electronic structure calculations. 
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Why should the reahstic modehng work? It is not an idle question, because 
eventually we would like to start with the electronic structure in the local-density 
approximation, construct the model by relying on this electronic structure, and apply 
this model for the strongly correlated systems. In the other words, we start with 
something what was claimed to be wrong for the strongly correlated systems and try to 
find its refinements by continuing to stay within this disputable picture and relying on 
these disputable approximations. Nevertheless, there are several reasons to believe that 
such a construction is possible and propose a meaningful strategy for improving LDA by 
staying within LDA and, at the same time, bypassing its limitations and shortcomings. 

(i) By the construction, the Kohn-Sham theory is the one-electron theory [2l E]. 
Therefore, it should provide a good starting point for the one-electron part of 
the model, no matter whether it is supplemented with LDA or not. Moreover, 
all contributions to the Kohn-Sham Hamiltonian coming from the exchange and 
correlations are local in the local-density approximation. Therefore, they can be 
easily subtracted in the process of construction of the model Hamiltonian in order 
to avoid the double-counting problem. 

(ii) The Coulomb correlations is certainly a weak point of LDA. However, LDA is an 
approximation in the density-functional theory, which is formulated for the ground 
state. Hence, it should be able to provide a good estimate for the total energy of 
the system [H [2], [3] . Then, there is always a chance to derive the effective Coulomb 
interactions from the total energy difference (or any other quantity uniquely related 
with the total energy) by applying certain constraint conditions, which would 
simulate the redistribution of certain order parameters (typically, the electron 
density of the system) . 

The goal of this review article is to outline the main ideas, principles, and methods 
of the realistic modeling as well as to illustrate abilities and perspectives of this direction 
for the solution of several practical questions related with the understanding of the real 
physical properties of transition-metal oxides and other strongly correlated systems. 
After specifying the model in Section Ej in the next two Sections we will discuss 
how the parameters of this model can be derived from the first-principles electronic 
structure calculations. Section |3] will be devoted to the construction of the one-electron 
Hamiltonian starting from the LDA band structure, while Section H] will deal with the 
problems of screening and calculation of the effective Coulomb interaction. We will 
try not only to illustrate the method, but also to present a clear physical picture for 
underlying ideas and approximation. Then, Section |5] will briefly summarize the methods 
of solution of the model Hamiltonian. In Section [6] we will give practical examples and 
show applications for realistic materials. We will derive the parameters of the model 
Hamiltonian and discuss what these parameters mean for understanding the electronic 
and magnetic properties of the considered compounds. Section [7] will contain brief 
summary and concluding remarks. 
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2. Model Hamiltonian 

First of all, we would like to specify our model for the the low-energy bands. We will 
be mainly dealing with the effective multi-orbital Hubbard model, 

"^ — 2^ 2-^ "'RR''^Ra'^R'/3 + 9 2-^ 2-^ ^a/37'5'^Ra'^R7'^R/3'^R<5' \^) 

RR' a/3 R afSjS 

where c{^„ (crq,) creates (annihilates) an electron in the Wannier orbital W^ of (typically, 
the transition-metal) center R, and a is a joint index, incorporating all remaining 
(spin and orbital) degrees of freedom, unless it is specified otherwise. The one- 
electron Hamiltonian /iRR' = ||/i^'j^/|| usually includes the following contributions: the 
site-diagonal part (R=R') describes the local level-splitting, caused by the crystal field 
and/or the relativistic spin-orbit interaction, whereas the off-diagonal part (R^^R') 
stands for transfer integrals, describing the kinetic energy of electrons. 

Uo.^,s = I dr I dr'W^^{r)Wl{r)v,Ur,r')Wl\r')Wi{r') ^ {WlW^\v,,,\Wlwi) 

are the matrix elements of the screened Coulomb interaction fscr(r,r'), which are 
supposed to be diagonal with respect to the site indices {R}. In principle, Ua/S'ys 
can also depend on the site-index R. Nevertheless, for the sake of simplicity of our 
notations, here and throughout in this paper we drop the index R in the notation of the 
Coulomb matrix elements. The intersite matrix elements of Uais-ys are typically small in 
comparison with the on-site ones. 

3. One-Electron Hamiltonian and Wannier functions 

The one-electron part of the model Hamiltonian (|3]) is typically identified with the Kohn- 
Sham Hamiltonian in the basis of Wannier functions representing the low-energy part 
of the spectrum [12]. Therefore, the concept and definition of the Wannier functions 
is one of the key parts of the methods, and we would like to start our discussion by 
making several general comments about the relationship between Wannier functions and 
localized atomic orbitals, which represent the basis of many computational schemes. 

Let us assume that there is a certain set of localized orbitals {xr} centered at 
the atomic sites {R} and specified by the orbital indices {a}a The corresponding 
Kohn-Sham Hamiltonian in the basis of {xr} will be denoted as H. The orbitals 
are orthonormalized and form a complete basis in the valence part of the spectrum, 
so that each eigenvector tpi of H can be expressed as a linear combination of {Xr}- 
The concrete examples of such bases can be the orthonormalized atomic orbitals or the 
muffin-tin orbitals ITSllTilfTS]. 



^ In the following, a set of nonorthonornialized atomic-like orbitals will be denoted as {Xr}- 
The orthonormalized orbitals, constructed from {xr,}i ai'e denoted as {xr}- Generally, such an 
orthonormalization can be performed numerically. 
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Since the Wannier functions {W^} are also defined as certain set of localized orbitals 
representing {ipi} [121 US], we immediately recognize that for the full Hamiltonian H, 
{Xr} can be regarded as one of the possible (and fully legitimate) choices for {Wr}. This 
is a natural result and advantage of the basis of localized atomic orbitals. In the plane- 
wave basis, the localized Wannier functions can be constructed from the eigenstates of 
H in the valence part of the spectrum, for example, by minimizing the square of the 
position operator (r^) [16]. However, we would like emphasize that this is nothing but 
an elegant way of constructing the compact atomic-like orbitals from the extended plane 
waves, a step which becomes rather unnecessary if one works from the very beginning 
in the atomic basis. 

However, what we typically need in the process of construction of the model 
Hamiltonians is different. For example, the solution of the many-body problem is 
practically impossible in the Hilbert space of states {xr.} of the full Hamiltonian H. 
Instead, one would like to concentrate on the behavior of a small number of L-bands, 
typically located near the Fermi level, and construct the Wannier basis only for this 
group of bands, which would be also orthogonal to the rest of the eigenstates of H. 
This causes an additional complication because the basis functions {xr}; though can 
be regarded as the Wannier functions for the full Hamiltonian H, are no longer those 
for any subspace of H. 

At present, there are two methods, which are typically used to circumvent this 
problem and construct the Wannier functions for the subspace of H: the projector- 
operator method [Ml [HI [IHl [19] and the downfolding method [5l [201 [2T1 [22] . 

3.1. The Projector- Operator Method 

In the projector-operator method, each (nonorthonormalized) Wannier function is 
generated by projecting a trial basis function Ixr), centered at the site R, onto the 
L-bands: 



Wk) = P\rn), (4) 



where 



P = J2\A)m (5) 

is the projector-operator onto the L-bands, tpi is the eigenstate of H, and i is a joint index 
combining the band index and the position of the momentum k in the first Brillouin 
zone. The functions {IV;a} can be numerically orthonormalized, 

|W^A) = EK)['^-'^']'r'r, (6) 

R't' 

where '5=||iS^*j^|| is the overlap matrix, 

S'^'n = {Wi.\Wi)^{xUP\xk)- (7) 
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Then, the one-electron part of the model Hamiltonian (^ is defined by the matrix 
elements of H in the basis of these orthonormalized Wannier orbitals: 

<R' = {Wi\H\Wi,). (8) 

3.2. The Downfolding Method 

The conventional downfolding method also implies that the atomic basis can be divided 
into two parts: {xr}={Xr}©{Xr}) so that the low-energy part of the spectrum is 
mainly represented by the {xR,}-states, while {xr} is the rest of the basis states, 
which mainly contribute to the higher-energy part. Then, each eigenstate of H can 
be identically presented as the sum \ipi) = \ipj) + \ip^), where \ipj) and \ipi) are expanded 
over the basis states of the "t" and "r" types, respectively. In this case, the Schrodinger 
equation for {ipi) takes the following form: 

(H'' - coM) + H'W,) =0 (9) 

H^'l^l) + {R'"- - uM) = 0, (10) 

where if *(''")*(^) are the blocks of matrix elements of H in the basis of "t" ( "r" )-states. 
The effective cu-dependent Hamiltonian H^^g is obtained by expressing Jt/'J) from fITOl) . 

m = -{H^'^~ur'H^w, (11) 

and substituting into ([9]). This yields 

H^siuj) = {H'' -u)- H''{H"- - ujy^H'\ (12) 

which formally acts only on \ipj). However, lipj) is only a part of the eigenvector, which 
is not orthonormalized. Therefore, Hes{^) should be additionally transformed to an 
orthonormal representation: 

h{uj) = S-^/^{uj)H,s{uj)S-'/\uj) + 00, (13) 

which is specified by the overlap matrix, 

S{uj) = 1 + H*''(H"- - ujY^H'K (14) 

The latter is obtained after the substitution of (ITTi) into the normalization condition: 

In the conventional downfolding method, h is typically evaluated in the center 
of gravity of the L-bands, ojq. Although the downfolding method does not explicitly 
require the construction of the Wannier functions, they can be formally reconstructed 
from hiuo) |22j . 

3.3. Downfolding as the Projector- Operator Method 

The conventional downfolding method is exact. However, this property is guaranteed by 
the cj-dependence of /i, which is hardly useful from the practical point of view. Formally, 
for each t/'j, w in (fT3|) should coincide with the eigenvalue of H corresponding to this ■j/'j. 
Moreover, h{u) retains an excessive information about if, so that the full spectrum of H 
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can be formally derived from h{uj). However, typically we do not need such a redundant 
information and would like to use h only for a small group of electronic states located 
near the Fermi level, and do it in the most exact form. 

For these purposes, the downfolding method can be reformulated as a projector- 
operator method and reduced to it [2^. The trick is to replace the original Hamiltonian 
H in the downfolding method by a modified Hamiltonian H', which has the same set 
of eigenvalues {Si} and eigenfunctions {ipi} in the region of L-bandso The rest of the 
eigenstates is not important for the construction of the one-electron part of the model 
Hamiltonian and can be placed to the region of infinite energies. Hence, we define H' 
in the following form: 

H' = Y, \i^i)e^m + eP± = PHP + eP^, (15) 

ieL 

where P±=l—P is the projector operator to the subspace orthogonal to the L-bands and 
e— j-oo. According to the choice of the basis functions {xr} and {xb} in the downfolding 
method, the latter mainly contribute to the high-energy part of the spectrum. Therefore, 
the overlap between i/ji in the low-energy part and any of {xr} should be small, so that 
all eigenvalues of [H'y^ are of the order of e. Then, it is intuitively clear that in the 
limit e— J-oo, the w-dependence in (ITTl) will disappear and the method will be exact. This 
statement can be proven rigorously, and such a proof can be found in [23] . 

Thus, by introducing the modified Hamiltonian of the form (IT5l) . the downfolding 
method can be naturally reformulated as the projector-operator method. The advantage 
of this procedure is that it allows us to go directly to the construction of the one-electron 
Hamiltonian h, and formally skip the step of the construction of the Wannier functions. 

Finally, we would like to note that ( !T5|) is nothing but a scissor-operator-like 
transformation of the original Hamiltonian H, and this strategy is different from the 
order-A^ muffin-tin orbital method, which was also used for the construction of the 
Wannier functions and where the basic idea was to make some approximations for the 
cj- dependence of the downfolded Hamiltonian [21] . 



3.4- Choice of Trial Orbitals and Localization of the Wannier Functions 

At the beginning of this section we have argued that the basis functions {xr} can be 
regarded as the Wannier functions of the full Hamiltonian H. Now, let us assume that 
each basis function is localized around the central atomic site and satisfies certain criteria 
of the "maximal localization", such that any linear combination of {xb} ^^^^ t)e "less 
localized" in comparison with the basis function Xn from the original basis set, or at least 
has the same degree of the localization, if we are dealing with the linear combination of 
orbitals centered at the same atomic sitejj However, this is not necessarily true if one 

•^ This procedure was already used in [5l [20l [21] [22] . However, the details have been explained only 
in [23]. 

The precise criterion of the maximal localization is not really important at this stage, because this 
is merely a mathematical construction and depending on the considered physical property one can 
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wants to construct the Wannier functions only for some part of the electronic structure, 
which is specified by the subspace L of the full Hamiltonian H. Due to the additional 
orthogonality condition to other bands, such a Wannier function will inevitably be a 
linear combination of {Xr}- Therefore, it will inevitably be less localized in comparison 
with the trial function Xb.- 

Nevertheless, one can try to "minimize" the deviation from the original formulation 
for the full Hamiltonian H and ask which single atomic orbital centered at the site R 
will be the best representation for the Wannier orbital. Therefore, we search a new set 
of trial functions in the form: 

a 

and find the coefficients {c^}, which maximize the projection (</'r|W^r[</'r,]) of 0^ onto 
the nonorthonormalized Wannier function constructed from 0^ using the projector- 
operator technique, |W^r[0r])=-P|0r)- It will automatically guarantee that 0^ is the 
best single-orbital representation for W^ in the projector-operator method among the 
trial orbitals of the form (IT6ll . By substituting W^r[0r] into the projection (0r|W^r[0r]), 
problem is reduced to the maximization of 



D = max{(4|P|4) - A((4|0*i,) - 1)| 



with respect to {c^}, where the Lagrange multipliers {A} enforce the orthonormality 
condition for {0r}. Then, the maximization of D is equivalent to the diagonalization 
of -Prr=||(Xr|-P|Xr)||) which is nothing but the site-diagonal part of the density 
matrix constructed from the L-bands in the basis of atomic orbitals {xr}- After the 
diagonahzation, we should simply pick up n eigenstates {0r}, corresponding to maximal 
eigenvalues {A}, where n is the number of Wannier functions centered at the atomic site 
HU These {0^} will maximize D. This procedure has been proposed in [20] without 
any proof. Then, some intuitive arguments have been given in [22]. The rigorous proof 
has been given in 



3. 5. Other parameters of model Hamiltonian 

Electronic structure calculations are typically supplemented with some additional 
approximations, like the atomic-sphere-approximation (ASA), where the LDA potential 

introduce different criteria of tire "maximal localization" . 

^ This can be paraphrased in a different way [5 . Of course, any set of the Wannier functions should 
be able to reproduce the total density at the site R. Each Wannier function consists of the central part 
(or "head"), located at the site R, and the tail, spreading over the other sites. By identifying {0^} 
with eigenstates corresponding to the maximal eigenvalues of the density matrix {A}, we guarantee 
that the main part of the density at the site R is described by the "heads" of the Wannier functions. 
Then, the remaining part of the density, corresponding to other (small) {A}, is described by the tails 
of the Wannier functions coming from other sites. This implies that the weights of these tails should 
be also small. 
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is spherically averaged inside certain atomic spheres [13] • If ASA is used in the process 
of construction of the model Hamiltonian, the parameters {h'^J^^,} will include all effects 
of the covalent mixing effects between atomic orbitals. However, there will be other 
effects, which are not accounted by {h^^,}. The most important contribution comes 
from the nonsphericity (n-s) of the Madelung potential for the electron-ion interactions, 
which contributes to the crystal-field splitting [25] |j The proper correction to ||/iRpi/|| 
can be computed in the Wannier basis as: 

^""-'h^A = E <^^' |b1?-bJ ^h). (17) 

RYR l^ + i ^1 

where Z^, is the total charge associated with the site R' (namely, the nuclear charge 
minus the screening electronic charge encircled by the atomic sphere), and r is the 
position of electron in the sphere R. 

The main idea behind this treatment is based on certain hierarchy of interactions in 
solids. It implies that the strongest interaction, which leads to the energetic separation 
of the L-band from other bands (Figure [1]), is due to the covalent mixing. For example, 
in many transition-metal oxides this interaction is responsible for the famous splitting 
between transition-metal t2g and Cg bands [26]. The nonsphericity of the Madelung 
potential is considerably weaker than this splitting. However, it can be comparable 
with the covalent mixing in the narrow L-band. Therefore, the basic idea is to treat 
this nonsphericity as a pseudo-perturbation [T3|, and calculate the matrix elements of 
the Madelung potential in the basis of Wannier functions constructed for spherically 
averaged ASA potential. 

The same strategy can be applied to the spin-orbit (s-o) interaction, which yields 
the following correction to ||/iR^||: 

Here, V is the self-consistent LDA potential and a is the vector of Pauli matrices. 

4. Effective Coulomb Interactions 

Generally, the matrix elements of the effective Coulomb interaction in the L-band are 
defined as the energy cost for moving an electron from one Wannier orbital, say W-^,, 

^ In fact, the contribution of Coulomb interactions to the crystal-field splitting is a tricky issue. Despite 
an apparent simplicity of the problem, one should clearly distinguish different contributions and not 
to include them twice, for example, to the one-electron and Coulomb interaction parts of the model 
Hamiltonian (j3]). In this sense, the use of full-potential techniques does not automatically guarantee 
the right answer. Note that the nonsphericity of on-site interactions is explicitly included into the 
second part of the model Hamiltonian ([3]). Therefore, in order to not to include it twice, one should 
subtract the corresponding contributions to the one-electron part originating from the Coulomb and 
exchange-correlation potentials in LDA. 
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populated by n^/^ electrons, to another orbital, say W^, which was initially populated 
by n-^a electrons |27] : 

Uaa/313 = E [tIr^ + 1, nji>p - l] - E [n^a, ?^R'/3] • (18) 

For Rt^R', the above matrix elements define the on-site Coulomb interactions, which 
are screened by intersite interactions. In principle, by considering different combinations 
of R and R', one can extract individual parameters of on-site and intersite interactions 
[22] . However, in practice, the intersite interactions are substantially smaller than 
the on-site ones, and can be neglected in many caseslH The total energy difference 
(TT8|) corresponding to R=R', but with different orbital indices a and /3, define the 
nonspherical part of on-site interactions, which is responsible for Hund's rules. More 
generally, one can consider an unitary transformation of a and (3 at each site of the 
system to a new basis, and repeat the same arguments. In the original basis, this would 
correspond to the transfer of an electron from a linear combination of Wannier orbitals 
at the site R' to a linear combination of Wannier orbitals at the site R. This procedure 
will define the full matrix of screened Coulomb interactions ?7=||f/a/375||. 

4.I. Constrained density-functional theory 

The dependence of the total energy £'[{^^0,}] on the individual occupation numbers 
{n-Ra} can be obtained by applying the constrained density-functional theory [29l [30l 



E [{nua}] = E[p] + J2 Vna {{W^\p\W^) - nR„) , (19) 

Ra 

where p(r, r') = Xli /«V^j (^)V^i('^') i^ ^^^ density operator constructed from the Kohn- 
Sham orbitals, p(r) = p(r, r) is the electron density, and {Vrq,} are the external 
potentials, which play the role of Lagrange multipliers and enforce the occupations 
of the Wannier orbitals to be equal to {wrq,}. As it was already pointed out in the 
Introduction, in spite of many limitations for the strongly correlated systems, LDA is 
formulated as an approximation to the theory of the ground-state. Therefore, there is 
always a hope that it will provide a good estimate for Uaapp as long as the latter is 
expressed in terms of the total energy difference flTSj) or any other quantity, which is 
uniquely related with this total energy difference. 

The total energy difference (IT8l) is typically replaced by the difference of Kohn- 
Sham eigenvalues {eRa} calculated for an intermediate configuration, by using Slater's 

Yet, one exception can be the parameters derived for the single-orbital model. In this case, the 
number of variables available for the localization of the Wannier functions is limited so that the latter can 
be rather extended in the real space [551. Therefore, the hare intersite Coulomb and exchange integrals, 
calculated in the basis of these Wannier functions, are no longer negligible [T71HH]- Nevertheless, at 
present it is not entirely clear to which extent these interactions in the L-band will be screened by other 
bands of the system. 
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transition state argumentslfl 

It implies the validity of Janak's theorem [33] 

ena = OE/nna- (21) 

However, in order to use this theorem, {W^r} must be the eigenvectors of Kohn-Sham 
equations ([1]). Otherwise, Janak's theorem is not applicable, and the eigenvalues {^Ra} 
are ill-definedjj Of course, this assumption does not work for the regular Wannier 
functions {VTr} constructed for the LDA Hamiltonian H\^^\ 

Nevertheless, there is a transparent and physical solution of this problem, which 
rehabilitates the use of Janak's theorem. Indeed, since the kinetic-energy term is 
explicitly included into the Hubbard model ([3]), it should not contribute to the total 
energy difference ( 1T8|) . Otherwise, it would be included twice, and we would face the 
problem of the double counting [30l [3T] . This can be done by artificially switching off 
all matrix elements of H between W^ and W^/, which is equivalent to switching off 
the transfer integrals. Another solution is to modify the Wannier functions {IVr} — )■ 
{^r} (apparently by considering another boundary conditions) in order to enforce the 
property {W^\H\W^,)=0 for the given LDA Hamiltonian H [22]0 These {W^} can be 
regarded as the eigenf unctions of H, that justifies the use of Janak's theorem. 

Finally, in the first order of ±|, ( l20ll can be transformed to 

Uaal3l3 = ds-Ra/ du-R^a, (22) 

where the derivative is calculated under the condition that the total number of electrons 
is conserved: nRa+^R'/3=const. 

It is important that in the process of solution of the Kohn-Sham equations, the 
Wannier orbitals (and all remaining electronic structure) are allowed to relax in order to 
to adjust the change of the occupation numbers {nua}- This relaxation is an important 
channel of screening. As we shall see below, the values of the screened Coulomb 

^ Here, we also use the reciprocity condition eu.'i3[nRa+^,n^'i3 — ^] = £rq["^Rq — ^,f^R'/3 + ^] for two 
fluctuations around the uniform state with riRa = nR'i3- 

^ For example, by defining erq as the diagonal matrix element of the Kohn-Sham Hamiltonian 
H, SRa — iW-^lHlW-^), does not necessary guarantee that this £rq is equal to dE/n-Ra, and (f20|) 
is consistent with the more general definition (J18p . 

^" Note that if {M^^} were the eigenvectors of H, all transfer integrals, which are defined as off-diagonal 
matrix elements of H with respect to W^ and W-^, would vanish (see Section I5TT|) . 
^^ This can be naturally done by reconstructing the Wannier functions from the matrix elements h'^^, 
derived from the downfolded method [22 . Nevertheless, it seems that for many strongly correlated 
systems, {W|£} do not substantially differ from {W|^}. The intuitive reason for that is that the 
transfer integrals in the narrow L-band are typically much smaller in comparison with other effects of 
the covalent mixing, which lead to the separation of the L-band from other bands and mainly determine 
the shape of the Wannier functions. Another reason is that for strongly correlated systems, the transfer 
integrals are typically much smaller than the on-site Coulomb interactions. Therefore, to certain extent 
it is unimportant whether the small transfer integrals are included in the definition of the large Coulomb 
interactions or not as they cause only small change of these interactions. 
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interactions in the L-band are defined not only by the extension of the Wannier functions 
in the ground state, but mainly by the change of these functions during the reaction 
{n-Ra,n-R'i3) ^ {nua+^yn-R'is—l). Thus, in order to calculate U, it is not sufficient to 
know the Wannier functions in the ground state. Even more important question is how 
these Wannier functions will be modified in the process of transfer of an electron from 
one site of the system to another. 

4-2. Random-Phase Approximation 

Another way of calculating the screened Coulomb interaction is to use the random-phase 
approximation (RPA), which constitutes the basis of the GW method [311 ESI |37]. In 
this case, the effective Coulomb interaction is calculated in terms of the response onto 
the external perturbation of the electron density, (Jpextlr), which is introduced as a probe 
of the electron system. RPA has many similarities with the constrained DFT. It consists 
of the following steps. 

(i) 5pcxt(i") creates the Coulomb potential ^Kxtlr) = e^ / (ir'|r— r'|^^(5pcxt(i"')5 which 
is similar to {Vaa} in the constrained DFT, and 5pext(r) itself has a meaning of 
the electron density, which is controlled by the orbital occupations {riRajo In 
principle, the perturbation of the electron density can also depend on time, that 
would correspond to the time- dependent DFT [361 EZ] • Although such an extension 
can be certainly done, we will consider only the static (or time-independent) limit, 
strictly following the definition ( lT8l) . 

(ii) 5\4xt(r) is treated as a perturbation. Then, the updates for the Kohn-Sham orbitals, 
{Sipi}, can be easily calculated in the first order of (5V^xt(r), by using the regular 
perturbation theory. The screening of the Coulomb interaction in RPA is entirely 
related with the change of {ipi} or, in the other words, the relaxation of the Kohn- 
Sham orbitals. 

(iii) The change of the electron density is calculated from {Sipi} in the first order of 
SVe^tir), namely, 6p{r) = Ei[^V'l(r)V^i(r) + V'|(r)(5V'i(r)]. 

(iv) Sp{r) results in the additional change of the Coulomb potential, also in the first 
order of 6Vcx_t{r): 5^H(r) = e^ / (ir'|r— r'|~^(5p(r'). 

(v) The steps (ii)-(iv) are repeated for the new potential, SVext{'i')+SV}i{r). Then, the 
potential is updated again, and so on. This is equivalent to the iterative solution 
of the Kohn-Sham equations within DFT. The constrained DFT also takes into 
account the change of the exchange-correlation potential, ^xc(r). However, this 
change is typically treated in LDA, where the effect is small. 

^^Yet, there is also a difference: JT4xt(r) is the local potential, whereas {14iq} act on the individual 
Wannier functions {M^^}. This corresponds to the external potential SVcxt{^, r') = J^Ra ^R-al^R)(^Rl 
in the form of the projector operators, which are essentially nonlocal. This is the one of the basic 
differences between constrained DFT and RPA. The role of this nonlocality in the screening of the 
Coulomb interactions is not fully understood 1381. 
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(vi) After summing up an infinite number of updates for 6Vii{r), the screened Coulomb 
interaction in RPA is defined as f/(r, r') = 5V^xt(r)/'^Pcxt(r') + ^^(i")/^Pcxt(i"')- 

Tlien, one can easily to show that f/(r, r') satisfies the Dyson equation [Ml 135] : 

[/(r,r') =e>-r'|-i+ f dr" j dv"'e^\v - T"\-^V{v\r"')U{v'\v'), (23) 

where 'P(r, r') is the polarization function, which is obtained from the first-order 
perturbation theory expansion for {5ipi}: 

The w- dependence of V{r,r',ijj) corresponds to the time- dependent perturbation 
(5t>ext(i") exp(— icjt). The static limit corresponds to P(r, r',0) = P(r, r'). The prefactor 
"2" stands for two spin channels in the case of non-spin-polarized LDA. 

4-3. Combining Constrained DFT and RPA 

The screened Coulomb interactions in RPA satisfies one very important property |39j . 
which directly follows from the Dyson equation (123!) . Suppose that one can identify two 
different channels of screening, by dividing the full polarization function (IMI) into two 
parts: P = Pi + P2- Then, the screened Coulomb interaction U can be computed in 
two steps. 

(i) By replacing total V by Vi in (125]) . one can compute the effective interaction Ui, 
which takes into account only the first channel of screening. 

(ii) Then, the final interaction U, corresponding to V1+V21 can be obtained from Ui 
again by using the Dyson equation (l23!l . but where the full polarization function is 
replaced by V2'- 

U{r,r') = t/i(r,r') + f dr" f dr"'Ui{r,r")V2{r",r"')U{r"' ,r'). 

This property has two important consequences. 

(i) In order to calculate U, we would like to start with the LDA band structure. 
However, for Mott insulators, LDA yields an incorrect metallic behavior in the 
region of L-bands. Therefore, we should get rid of this unphysical metallic screening, 
which will always appear in RPA if one starts from LDA. More generally, RPA may 
not be a good approximation for treating the electron correlations in the narrow 
L-bands. Therefore, the general strategy is to get rid of the RPA self-screening 
associated with all kind of transitions between L-bands in the polarization function 
(plj) . and to replace it by a more rigorous model of screening in the process of 
solution of the Hubbard model ([3]) [211 EH]- The suppression of the transitions 
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between L-bands in RPA is similar to switching off all transfer integrals between 
Wannier functions in the constrained DFT [22]J^n 



(ii) Our goal is to evaluate the screening of the Coulomb interactions associated with 
the relaxation of the Wannier functions {W^r} in the region of L-bands as well as of 
all other Kohn-Sham orbitals {tpi} in the remaining part of the spectrum. Suppose 
that we are working in a (flexible) basis of atomic orbitals {xr}; and both {M^r} 
and {ipi} are expanded over this basis: 

?}-5:4.xi.. (25, 

^* J R'/3 

Then, within linear response theory, the relaxation of {VF^} ^^^ {V'j} "wiH consist 
of two parts: the relaxation of the basis functions {xb} ^^^ ^^^ relaxation of the 
coefficients {c^} of the expansion over these basis functions, which corresponds to 
the change of the hybridization in the process of screening. Therefore, it has been 
proposed (though somewhat heuristically) that the screened Coulomb interactions 
in the L-band can be computed in two steps [221 EH] • For the sake of clarity, let 
us assume the main contribution to the L-bands comes from the 3d orbitals of the 
transition-metal sites. 

(a) The first step takes into account the screening caused by the relaxation 
of the basis functions {Xr}- It can be easily evaluated in the framework 
of constrained DFT. In this case, the Kohn-Sham equations are solved 
numerically on a finite grid of points in the real space, and {n-Ra} are the 
constrained 3d occupation numbers inside certain atomic spheres [29]. On 
the other hand, since RPA is based on the perturbation theory expansion, 
similar calculations in terms of the polarization function (^^ would require an 
enormous number of states {ipi}, both in the occupied and unoccupied parts of 
the spectrum. Probably, this is one of the reasons why the conventional RPA 
calculations converge very slowly with the increase of the number of bands in 
the unoccupied part of the spectrum [IQ], and are extremely sensitive to the 
treatment of the core states [35j. The screening caused by the relaxation of 
the Kohn-Sham orbitals {ipi} in the subspace orthogonal to L (for example, 
the screening of localized 3d electrons by outer electrons in the transition- 
metal compounds) can be also easily taken into account in the conventional 
constrained DFT calculations [30l |3ll [32] . 

(b) However, what the constrained DFT cannot do is to treat the screening of 
Coulomb interactions amongst 3d electrons by the same 3d electrons [381 SI] o 

^^ Indeed, suppose that we have replaced {ipi} for the L-bands in (|24l) by the Wannier functions {W^}, 
and switched off all transfer integrals. In the ground state of such an auxiliary system, all occupation 
numbers are equal, ft—fj. Therefore, the polarization function (j24p will vanish. 

^^ The procedure implies that the atomic 3d charges can be divided into two parts, where the first 
part is subjected to the constraint conditions, while the second part is allowed to participate in the 
screening. Although, it can be formally done within constrained DFT [42], the actual computational 



Combining DFT and Many- Body Methods to Understand Correlated Materials 17 

Since the atomic 3d occupation numbers are rigidly fixed by the constraint 
conditions, the 3d electrons from other sites of the system cannot compensate 
the change of the number of the 3d electrons at the central site, and vise 
versa. However, such a channel of screening may exist. Suppose that 
our L-band is mainly constructed from the transition-metal 3d orbitals (the 
concrete example is the t2g band in the transition-metal oxides), and there is 
another, say oxygen 2p band, which has an appreciable weight of the atomic 
3d orbitals (of both t2g and Cg symmetry) coming from the hybridization 
between oxygen and transition- metal sites (Figure [2]). Furthermore, the 
redistribution of the electron density in the t2g band, associated with the 
reaction (?T.Ra,'n.R'^)^(riRa+l,?T.Ft/^— 1), will change the Coulomb potential 
around each transition-metal site. If the number of t2g electrons is increased by 
the constraint conditions, the Coulomb potential becomes more repulsive and 
vise versa. The more repulsive Coulomb potential will additionally push the 
3d states from the oxygen 2p band to the higher-energy part of the spectrum. 
Therefore, around certain transition-metal sites, the change of the number of 
3d electrons in the t2g band will be partly compensated by the change of the 
3(i-electron density in the region of oxygen 2p band. This channel of screening 
can be easily evaluated in RPA, by rewriting f l23|) in the matrix form 

U =\l- uv\ u, (26) 

and assuming that all other channels of screening are already included in the 
definition the "bare Coulomb interaction" u, derived from the constrained DFT 
|22j . Since the polarization matrix V in ( l26l) is aimed to describe the self- 
screening of the 3d electrons, it should consist of the matrix elements of ( 12^ in 
the basis of atomic 3d orbitals, after the subtraction of the unphysical metallic 
screening associated with RPA transitions between t2g band. 

5. Solution of Model Hamiltonian 

It is virtually impossible to provide a comprehensive analysis of all possible methods of 
the solution of the low-energy model ([3]), and this is definitely beyond the scopes of this 
review article. 

One option is the dynamical mean- field theory (DMFT) |13], which becomes one 
of the popular low-energy solvers today. The idea of DMFT is to map the many- 
body lattice problem to a single-site impurity problem with effective parameters. The 
vast majority of DMFT applications for realistic compounds have been focusing on 
the analysis of spectroscopic properties, especially in the context of the metal-insulator 
transition, although some extensions for the ground-state properties, such as calculations 
of the total energies and phonons, are also available today. Many examples of recent 

scheme is rather laborious, and in many cases the procedure of dividing the atomic charges into the 
"screened" and "screening" parts is not weh defined. 
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Density of States (arb. units) 

Figure 2. A schematic view on the change of the p-d hybridization in the oxygen 
2p and t2g bands of the transition-metal oxides associated with the repopulation of 
the Wannier orbitals at the neighboring transition-metal (TM) sites at two sides of the 
reaction 2((i")^(i"+^+(i"^^. Left panel corresponds to the ground-state configuration, 
2d". In the right panel, the removal (addition) of an electron from (to) the Wannier 
orbital in the t2g part of the spectrum is simulated by the shift of these orbitals relative 
to the Fermi level (shown by dot-dashed line). Around each transition- metal site, it 
changes the Coulomb potential, which controls the distribution of the Sd-states as well 
as the hybridization between transition-metal 3d and oxygen 2p states. Generally, the 
removal of an electron from (or the addition of an electron to) the Wannier orbital is 
partially compensated by the change of the amount of the 3d-states (shown by shaded 
area), which is admixed into the oxygen 2p band. This transfer of the spectral weight 
works as an efficient channel of screening of the local Coulomb interactions in the 
transition- metal oxides [22] . 



applications of DMFT can be found in the review articles pH |l5l HH]. The conventional 
DMFT becomes exact in the limit of infinite coordination numbers or, equivalently, 
infinite dimensions, when all nonlocal correlations vanish. In order to treat these 
nonlocal correlations, it is essential to go beyond the single-site approximation. This 
is one of the challenging problems in DMFT, and the recent progresses along this line 
can be found in |17]. Another limitation of DMFT is that, in order to be exact, it 
is typically used in the combination with the Quantum Monte Carlo (QMC) method, 
which provides an exact solution for the impurity model. However, current applications 
of the QMC method are typically restricted by rather high temperatures, which are 
substantially higher than, for example, the magnetic transition temperatures in many 
strongly correlated compounds. From this point of view, the method does not appear to 
be sufficiently useful for studying the phenomena of spin and orbital ordering. Probably, 
some of these difficulties may be overcome by applying the projective QMC method |18]. 
Unlike DMFT, the path-integral renormalization group (PIRG) method is mainly 
oriented on the description of the ground-state properties of strongly correlated systems. 
The entire procedure includes the following steps [HI jSOj [51]: 
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(i) The numerical construction of truncated basis of C Slater determinants, which 
provides the best representation for the ground-state wavefunction; 

(ii) Calculation of the total energy and its variance in the obtained basis; 

(iii) Extrapolation of the obtained results to the full Hilbert space, which is achieved 
by a systematic increase of C 

The PIRG method has been recently applied as the low-energy solver for studying the 
correlation effects in the t2g bands of Sr2V04 [21j and YVO3 [52] . 

In the rest of this section we will discuss some details of the solution of the model 
Hamiltonian (|3]), which will be directly used for applications considered in Section |6l 
We start with the simplest Hartree-Fock method, which totally neglects the correlation 
effects. Then, we consider simple corrections to the Hartree-Fock approximation, which 
include some of these effects. One is the perturbation theory for the total energy, and 
the other one is the variational superexchange theory. 

All model calculations are performed in the basis of Wannier functions {VFr}, which 
have a finite weight at the central transition-metal site as well as the oxygen and other 
atomic sites located in its neighborhood. In order to calculate the local quantities, 
associated with the transition-metal atoms, such as spin and orbital magnetic moments 
or the distribution of the electron density, the Wannier functions are expanded over the 
original basis {Xr}- Then, all aforementioned quantities are calculated by integrating 
over appropriate regions of the real space surrounding the transition-metal sites, like 
the atomic spheres in the LMTO method jlSj [T^ [T5] . 

5.1. Hartree-Fock Approximation 

The Hartree-Fock method provides the simplest approximation to the many-body 
problem (|3]). In this case, the trial many-electron wavefunction is searched in the form 
of a single Slater determinant \S{ipk}), constructed from the one-electron orbitals {^Pk}- 
In this notation, fc is a collective index combining the momentum k of the first Brillouin 
zone, the band number, and the spin (s= f or l) of the particle. The one-electron 
orbitals are subjected to the variational principle and requested to minimize the total 
energy 

EnF = mm{S{ipk}\n\S{ipk}) 

Wk} 

for a given number of particles A/". This minimization is equivalent to the solution of 
Hartree-Fock equations for {ipk}'- 

k + v]\^k)=ek\^k), (27) 
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where /ik^||/i^ || is the one-electron part of the model Hamiltonian (^ in the reciprocal 
space, /if =Er' h^-R'e~"^'^^~'^'\ and V=||V„/3|| is the Hartree-Fock potential0 

Va/3 = 2_^ (UalS^S — UaS^/s) ^^5- (28) 

7(5 

Equation ( 1271) is solved self-consistently together with the equation 

occ 
A; 

for the density matrix n=||no/3|| in the basis of Wannier functions. After iterative 
solution of the Hartree-Fock equations, the total energy can be computed as 



occ _. 



2 

k a(3 

By knowing {sk} and {fk}, one can construct the one-electron (retarded) Green 
function, 

U - Ek + iS 



W'(w) = Yl - ^' 



k 

which can be used for many applications. For example, the interatomic magnetic 
interactions corresponding to infinitesimal rotations of spin magnetic moments near 
the equilibrium can be computed as [531 EU • 

^RR' = ^ImJ ' dioTiL {gl^,{io)AVgi,^{uj)Av} , (29) 

where Q^^,=^Tis{{l±d'z)QB.-R.'} is the projection of the Green function onto the 
majority (t) and minority (J,) spin states, AV=TTs{crzV} is the magnetic (spin) part of 
the Hartree-Fock potential, Trg (Tr^) denotes the trace over the spin (orbital) indices, 
1 and az is the unity and Pauli matrix, respectively, and ep is the Fermi energy^ 

The parameters {Jrr'} are not universal, and depend on the magnetic state in 
which they are calculated, for example, through the change of the orbital ordering [5] 
or the change of the electronic structure by the magnetic ordering [71 [55] . 

^^ For the sake of simplicity, we drop the atomic index R in the notations of Vq^, although such 
a dependence can take place (for example, in the case of inequivalent transition-metal sites in the 
distorted perovskite structure) , and was actually taken into account in realistic calculations considered 
in Section [B] 

^^ According to the definition ([2^ . Jrr'>0 (<0) means that for a given magnetic state, the spin 
arrangement in the bond (RR ) corresponds to the local minimum (maximum) of the total energy. 
However, in the following wc will use the universal notations, according to which Jrr'>0 and <0 will 
stand the ferromagnetic and antifcrromagnetic coupling, respectively. 
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5.2. Second Order Perturbation Theory for the Correlation Energy 

The simplest way of going beyond the Hartree-Fock approximation is to include the 
correlation interactions in the second order of perturbation theory for the total energy 
[56| [571 EH]- It shares common problems of the regular (nondegenerate) perturbation 
theory. Nevertheless, by using this technique one can calculate relatively easily the 
corrections to the total energy, starting from the Hartree-Fock wavef unctions. This 
method is expected to work well for the systems where the orbital degeneracy is lifted (for 
example, by the crystal-field splitting) and the ground state is described reasonably well 
by a single Slater determinant, so that other corrections can be treated as a perturbation. 
The correlation interaction (or the interaction of fluctuations) is defined as the 
difference between true many-body Hamiltonian ([3]), and its one-electron counterpart, 
obtained at the level of the Hartree-Fock approximation: 

^^ = 5Z I 2 5Z f^"/3754a47CR/3Cm - Yl ^"/^cLcr/j I • (30) 

R V alS-yS a/B J 

It is important to note that although some of the matrix elements UajB-ys can be large, 
they also contribute to the Hartree-Fock potentials Vap- Therefore, generally, one can 
expect some cancelation of contributions in the first and second parts of (!30!) . which 
formally extend the applicability of the perturbation theory even for relatively large 
Ua/s-ys. For example, in a number of cases such a strategy can be applied even for the 
bare Coulomb interactions in isolated atoms [59] . 

By treating He as a perturbation, the correlation energy can be easily estimated 
as [561 [53 [58]: 



^ " V EMS)-EMG) ' ^ ^ 

where \G) and \S) are the Slater determinants corresponding to the low-energy ground 
state (in the Hartree-Fock approximation), and the excited state, respectively. Due 
to the variational properties of the Hartree-Fock method, the only processes that may 
contribute to E^j are the two-particle excitations, for which each of IS") is obtained 
from \G) by replacing two one-electron orbitals, say ipki and ipk2, from the occupied 
part of the spectrum by two unoccupied orbitals, say ipk^ and (pk^ [59j . Hence, using the 
notations of Section [2], the matrix elements take the following form: 

{SindG) = ihhlvscrlkih) - (hhlvscAhh) ■ (32) 

By employing further the approximation of noninteracting quasiparticles, the 
denominator in ([3T!) can be replaced by the linear combination of Hartree-Fock 
eigenvalues: ii^HF(5')— -Ehf(G') ~ ek3+ek4,—eki—£k2 [HSl [571 ISHj- The matrix elements 
dMD satisfy the following condition: {S\nc\G)^^Y.n(^'''^'^^'~^'~^'^''^ (^ being the 
number of sites), provided that the effective Coulomb interactions are diagonal with 
respect to the site indices. In the second-order perturbation theory one can estimate 
relatively easily both on-site (R=0) and intersite (Rt^O) contributions to E^' . The R=0 
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term corresponds to the commonly used single-site approximation for the correlation 
interactions, which becomes exact in the limit of infinite spacial dimensions [43j. 

In principle, one can go beyond the second order perturbation theory and consider, 
for example, the single-site approximation for the T- matrix |60]. In this case, the 
expression for the energy of electron-electron interactions has the same form as in the 
Hartree-Fock method, but with U being replaced by the effective T-matrix, which takes 
into account the correlation effects. The method has been employed for the series of 
distorted transition-metal perovskite oxides [61], where the degeneracy of the Hartree- 
Fock ground state is lifted by the crystal field. It that case, the application of the 
T-matrix theory changed only some quantitative conclusions, whereas the main trends 
for the correlation energy were captured already by the second order perturbation theory. 

5.3. Atomic Multiplet Structure and Superexchange Interactions 

The variational superexchange theory takes into account the multiplet structure of 
the excited atomic states. By using this technique one can study the effect of the 
electron correlations on the spin and orbital ordering. However, it is limited by 
typical approximations made in the theory of superexchange interactions, which treat 
all transfer integrals as a perturbation. 

The superexchange interaction in the bond (RR') is basically the gain of the kinetic 
energy, which is acquired by an electron at the center R in the process of virtual 
hoppings into the subspace of unoccupied orbitals at the center R', and vice versa 
[TT1[T0]. Therefore, the energy gain caused by virtual hoppings in the bond (RR') can 
be found as [5l[ini[62]: 



7r.r/ — — i G 



' r -C/R'M ' 



a , (33) 



where G is the ground-state wavefunction of the lattice of isolated centersjlH each 
of which accommodates n electrons, E-r^im and |R'M) stand for the eigenvalues and 
eigenvectors of the excited (ra+1 electron) configurations of the center R', and ^r' is a 
projector operator, which enforces the Pauli principle and suppresses any hoppings into 
the subspace of occupied orbitals at the center R' [3] . 

The formulation is extremely simple for the n = 1 compounds, like YTiOs and 
LaTiOs. In this case, there is only one electron residing at each transition- metal site. 
This is essentially an one-electron problem, where each atomic state is described by 
certain one-electron orbital (/Jr and G is the single Slater determinant constructed from 
{v^r} belonging to different transition- metal sites |5]. A similar formulation can be 
performed for the hole spin-orbitals {ctR} of compounds where at each lattice center 
there is only one unbalanced hole. Such a situation holds for the alkali hyperoxides, 
which will be considered in Section 16.41 

^^ In the present context, the "lattice center" means either isolated atomic site or a molecule. An 
example of the molecular solid will be considered in Section 16.41 
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The total energy of the system in the superexchange approximation is obtained after 
summation over all bonds, which should be combined with the site-diagonal elements, 
incorporating the effects of the crystal-field splitting and the relativistic spin-orbit 
interaction: 

EsE = y^(v?R|/^RR|'/?R) + ^ TrR'. 
R (RR) 

Finally, the set of occupied orbitals {v^r} is obtained by minimizing E^e- This can be 
done by using, for instance, the steepest descent method. 



6. Examples and Applications for Realistic Compounds 

6.1. Cubic Perovskites: SrVO^ 

SrVOa is a rare example of perovskite compounds, which crystallizes in the ideal cubic 
structure. It attracted a considerable attention in the connection with the bandwidth 
control of the metal- insulator transition [H |63]. The region of interest is the t2g band, 
which is located near the Fermi level (Figure E]). 




-6-4-20246 
Energy (eV) 



Figure 3. Total and partial densities of states of SrVOa in the local-density 
approximation. The shaded area shows the contributions of the vanadium 3d states. 
Other symbols show positions of the main bands. The Fermi level is at zero energy. 



6.1.1. Transfer integrals and Wannier functions. For cubic compounds, the separation 
of the basis functions into {xt} and {Xr}, which is required in the downfolding method, 
is rather straightforward: three t2g orbitals centered at each vanadium site of SrVOs 
are taken as the {xt} orbitals, and the rest of the basis functions is associated with the 
{Xr} orbitals. 

The downfolding procedure is nearly perfect and well reproduces the behavior of 
three t2g bands (Figure H]). As expected for cubic compounds, the nearest-neighbor ddn- 
interactions mediated by the oxygen 2p orbitals are the strongest (Tabled]). For the 



Combining DFT and Many-Body Methods to Understand Correlated Materials 



24 



xy-orbitals, it operates in the x and y directionsHj However, there is also an appreciable 
(i(i5-interaction operating in the "forbidden" direction (for example, the direction z in 
the case of the xy orbitals). These interactions are mediated by the strontium Ad 
orbitals and strongly depend on the proximity of the latter to the Fermi level. The 
transfer integrals connecting different t2g orbitals are small and contribute only to the 
longer-range interactions separated by the vectors (a, a, 0) and {a, a, a), where a is the 
cubic lattice parameter [64j. Other interactions are considerably smaller. 
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Figure 4. Left panel: Crystal structure of cubic perovskites and atomic wavefunctions 
mediating transfer interactions between vanadium t2g orbitals. The standard Y{xy)- 
0(jj)-\(xy) and \{xy)-0{x)-\(xy) interactions operate in the x- and y-directions, 
respectively. The V(a;y)-Sr(3z^-r^)-V(a;j/) interaction operate in the "forbidden" z- 
direction. Right panel: LDA energy bands of SrVOa obtained in the original electronic- 
structure calculations using the LMTO method and after the "tight-binding" (TB) 
parametrization using the downfolding method. Notations of the high-symmetry points 
of the Brillouin zone are taken from 1651 . 



Table 1. Transfer integrals (in meV) between first, second, and third nearest neighbors 
in SrVOa, separated by the vectors (a, 0, 0), (a, a, 0), and (a, a, a), respectively (a being 
the cubic lattice parameter). 
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The shape of the Wannier functions is explained in Figure [Sl^^l Since t2g band is 
an antibonding combination of the atomic vanadium 3d-t2g and oxygen 2p orbitals, the 
Wannier function has nodes between vanadium and oxygen sites. Right panel of Figure 
|5] illustrates the spacial extension of the Wannier functions. It shows the weight of the 

^* Similar dependencies for the yz and zx orbitals are obtained by the cyclic permutation of the indices 
X, y, and z. 

^^ These Wannier functions have been reconstructed from the one-electron part of the downfolded 
Hamiltonian using the ideas of the LMTO method [13l [T4j [15] . The procedure has been explained in 
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Wannier function accumulated around the central vanadium site after adding every new 
sphere of the neighboring sites. Since Wannier functions are normalized, their total 
weight is equal to one. In the case of SrVOa, 77% of the this weight belongs to the 
central vanadium site, 16% is distributed over four neighboring oxygen sites, about 5% 
belongs to the next eight strontium sites, and 1% - to the eight oxygen sites located 
in the fourth coordination sphere. Other contributions are small. Another quantity, 
characterizing the spread of the Wannier functions, is the expectation value of square of 
the position operator, (r^) = (WR|(r— R)^|iyR,) |T6], which in the case of SrVOs is about 
2.37 A2 [22J0 
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Figure 5. Wannier function for tlie i2g band of SrVOs ^\. Left panel shows tlie 
contour plot of the ly-orbital in the (001) plane. The solid and dashed lines correspond 
to the positive and negative values of the Wannier function. The vanadium atom is 
located in the center of the plot, and other sites are four oxygen atoms located in 
its nearest neighborhood. Around each site, the Wannier function increases/decreses 
with the step 0.04 from the values indicated on the graph. Right panel shows the total 
weight of the Wannier function accumulated around the central vanadium site after 
adding every new sphere of the neighboring sites. 



For cubic perovskites, the transfer integrals can be extracted from first-principles 
electronic structure calculations in several different ways. For example, one can simply 
fit the LDA band structure in terms of the Slater-Koster parameters [6^ . However, the 
situation becomes increasingly complicated in materials with lower crystal symmetry, 
like distorted perovskite oxides, which will be considered below. First, the number of the 
Slater-Koster parameters, permitted by the symmetry, increases dramatically. Second, 
the form of these transfer integrals becomes more complicated and differs substantially 
from cubic compoundso Therefore, it seems that for complex systems the only way 

^'^ Somewhat smaller value (1.91 A^) has been reported in \Q6.. Some overestimation of (r^) is caused 
by some additional approximations used in the process of reconstruction of the Wannier functions from 
the downfolded one-electron Hamiltonian, which has been employed in |22j . Presumably, the direct 
application of the projector-operator method can do a better job. Nevertheless, as it was already 
pointed out in Section 13.31 after the transformation ()15p of the Kohn-Sham Hamiltonian, the transfer 
integrals derived from the downfolding method are totally equivalent to the ones obtained in the 
projector-operator method. Thus, although the Wannier functions reported in [22: may suffer from 
some additional approximations, the transfer integrals are essentially correct. 
^^ One example is the mixing of the t2g and eg orbitals by the orthorhombic distortion, which does not 
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to proceed is to use straightforward numerical algorithms, like the formal downfolding 
method. 

6.1.2. Effective Interactions. Applications of constrained DFT to the transition-metal 
oxides have been widely discussed in the literature [30l ETJ EHl IMl [70]. For example, the 
effective Coulomb interaction between Sd electrons in SrVOs can be computed in the 
following way [22] : 

(i) In the supercell geometry, one can introduce the "charge-density wave" , describing 
the modulation of the atomic 3(i occupations around the "ground-state" 
configuration with n = 1, n^ = n+6n cos{qR.) , where q is the propagation vector 
of the charge- density wave. 

(ii) Then, from the constrained DFT calculations, one can derive the Kohn-Sham 
eigenvalues {^r}, corresponding to this charge-density wave, and find the Fourier 
image of the effective Coulomb interaction as Uq = de^/dn^. 

(iii) Finally, the parameters of Coulomb interaction in the real space are obtained after 
the Fourier transformation of MqO 

For SrVOs, this procedure yields the following parameters of the on-site Coulomb 
interaction u = 10.1 eV and the nearest-neighbor Coulomb interaction v = 1.2 
eV. The intraatomic exchange interaction (j) can be derived by constraining the 3d 
magnetization density [29l [67] o This yields j = 1.0 eV. By knowing only u and j in 
the atomic limit, one can reconstruct the full 5x5x5x5 matrix ii of interactions between 
the 3d electrons, as it is typically done in the LDA+?7 method [71]. Some details of this 



procedure are explained in Appendix A 



In order to appreciate the magnitude of screening of different interaction parameters 
obtained in the constrained DFT, it is instructive to compare them with bare 
interactions. For example, the values of bare Coulomb and exchange integrals, calculated 
from 3d wavefunctions of the vanadium atoms, are 21.7 and 1.2 eV, respectively. The 
bare Coulomb interaction between neighboring vanadium sites, e^/a, is about 3.7 eV. 
Thus, in the constrained DFT, the on-site Coulomb interaction is reduced by factor 
two, the intersite Coulomb interaction is reduced by almost 70%, and the intra-atomic 
exchange interaction is reduced by 20%. All these interactions are further reduced by 
relaxation effects, related with the change of the hybridization. 

As it was already pointed out in Section 14. 3[ because of the hybridization, the 
transition-metal 3d states may have a significant weight in other bands. For example, in 
SrVOs besides the vanadium t2g band, the 3d states contribute to the vanadium eg as well 

occur in the cubic compounds. 

^^ For example, by considering only on-site (u) and nearest-neighbor intersite (w) interactions, we 
would have Uq — u — v J2r cos(qR), etc. 

^^ For example, if tor, is the 3d-magnetization, tor = n^—n^, and e^ and £^ are the Kohn- 
Sham eigenvalues for the majority- and minority-spin states, respectively, the parameter of intraatomic 
exchange interaction is given by j = 2d(e^—e^)/d'm-[i. 
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as to the oxygen 2p bands (Figure E])- If the number of t2g electrons changes, it causes 
some change of the Coulomb potential, which affects the distribution of the vanadium 
3d states in other parts of the spectrum. For example, if at certain vanadium site, 
the number of t2g electrons increases, the Coulomb potential becomes more repulsive. 
Therefore, the 3d states of this vanadium site will be pushed from the oxygen 2p band 
to a higher energy region. This causes some change of the coefficients {c^} of the 
expansion of the Kohn-Sham orbitals over the basis functions fl25l) or the change of the 
hybridization. This mechanism is responsible for an additional channel of screening of 
Coulomb interactions, which can be evaluated within RPA. In these calculations, the 
matrix ii, obtained in the constrained DFT method is used as the starting point, while 
the RPA itself is employed in order to evaluate the screening of 3d interactions in the 
vanadium t2g band by the same 3d states, which contribute to other bands. Thus, the 
problem is reduced to evaluation of the 3d matrix elements of the polarization function 
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Figure 6. Left panel shows the local density of states of SrVOa with the notation of 
the main interband transitions, which contribute to the polarization function in RPA: 
0{2p)^Y{eg) (1), 0{2p)^Y(t2g) (2), and V(i2ff)^V(eg) (3). Right panel shows the 
intraorbital Coulomb interaction U as obtained in the constrained DFT (denoted as 
"c") and after including the screening caused by different interband transitions in RPA. 
The RPA results show the screening corresponding to each type of transitions in the 
polarization function as well as the final value of U, which incorporates the effect of 
all three transitions [22]. 



According to the electronic structure of SrVOs, one can identify three main 
contributions to the polarization function, associated with the following interband 
transitions: oxygen 2p — )■ vanadium Cg, oxygen 2p — )■ vanadium t2g, and vanadium 
t2g — )■ vanadium Cg. 

The details of RPA screening are explained in Figure El For these purposes, it 
is convenient to introduce three Kanamori parameters [60]: the intraorbital Coulomb 
interaction 

U= I rfry"drXyr)W^.,(r)t;,„(r,r')W^4(r')W^.,(r'), 
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the interorbital Coulomb interaction 

W = j dv j dv'Wly{v)W,y{v)v,,,{ry)Wl{v')Wy,{v'), 
and the exchange interaction 

J = j dv j dv'Wl{v)Wy,{v)vUvy)Wly{v')Wy^{v'). 

In the atomic hmit, all interactions between t2g electrons are reduced to either lA, 
W, or J', and there is no other types of interactions connecting the t2g orbitals (see 



Appendix A ). Below we will argue that similar property holds even after the RPA 
screening. 

In addition to the final value of U, Figure [6] shows the screened interactions 
corresponding to each type of transitions in the polarization function. The screening 
caused by the change of the hybridization is very efficient. For example, in comparison 
with the constrained DFT, the intraorbital interaction U is reduced from 11.2 to 2.5 eV 
(i.e., by more than factor four). The main contribution to this screening comes from the 
oxygen 2p — )■ vanadium Cg and oxygen 2p — )■ vanadium t2g interband transitions in the 
polarization functions. Since the hybridization between vanadium t2g and Cg orbitals is 
small in perovskite compounds with the simple cubic structure, the screening associated 
with the transitions between vanadium t2g and Cg bands is also small. 

The dependence of the screened Coulomb interactions on the number of electrons, 
nt2gj accommodated in the t2g band is shown in Figured [22]. The calculations have 
been performed in the rigid-band approximation and using the electronic structure of 
SrVOs. Such an analysis may be useful for understanding the doping-dependence of the 
effective Coulomb interactions. 
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Figure 7. Doping-dependence of Kanamori parameters in SrVOs [22] : the intraorbital 
Coulomb interaction U, the interorbital Coulomb interaction U' , and the exchange 
interaction J^ versus the number of electrons in the t2g band, nt2g ■ 



The Coulomb interactions reveal a monotonic behavior as the function of doping. 
The screening is the most efficient when the whole t2g band is empty {ut^g = 0). 
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The situation corresponds to SrTiOs, where all transitions from the oxygen 2p to the 
transition- metal t2g band can contribute to the screening (see Figure |6]). This channel 
of screening is closed when the t2g band is filled [rit^g = 6). In the latter case, only the 
oxygen 2p — )■ transition-metal Cg interband transitions may contribute to the screening. 
Hence, the effective Coulomb interaction becomes large. 

The screening of the exchange integral J^ practically does not depend on the 
doping. The Kanamori rule, U = W + 1J ^ which was originally established for the 
spherical environment in isolated atoms, works well also for the t2g manifold in the 
cubic compounds, even after the screening of t2g interactions by other electrons. 

This result support an old empirical rule suggesting that only the Coulomb integral 
U is sensitive to the crystal environment in solids. The nonspherical interactions, which 
are also responsible for Hund's first and second rules, appears to be much closer to their 
atomic values and practically insensitive to the screening [72l [73l [71]. 

It is important to note that the obtained values of effective Coulomb interactions 
are substantially smaller than the experimental parameters derived from the analysis 
of photoemission spectra [751 ES]- However, this is to be expected. Note that the 
photoemission spectra are typically interpreted in the cluster model, which treats 
explicitly all transition-metal 'id as well as the oxygen 2p states. However, in the model 
([3]) we would like to keep only the transition-metal t^g bands and include the effect of 
other bands implicitly^ i.e. through the renormalization of interaction parameters in the 
t2g band. Therefore, our parameters should be generally smaller in comparison with the 
ones derived from the cluster model. As it was already discussed above, the transfer 
of an electron, associated with the reaction {n^a-,'n-^ j^) ^ ('^Ra+l,''^R'/3— 1) will cause 
some change of the electronic structure in the region of oxygen 2p and transition-metal 
Cg bands, which tends to compensate the change of the number of the 3(i electrons in the 
t2g band. Since the oxygen 2p and transition-metal Cg bands are eliminated in our t^g 
model, this change of the electronic structure is effectively included into the screening 
of Coulomb interactions in the t^g band, that naturally explains smaller values of the 
parameter U. 

Finally, the obtained value of intraorbital Coulomb interaction lA = 2.53 eV is 
substantially smaller than U ?» 5.55 eV, which is typically used in DMFT calculations 
in order to reproduce the experimental photoemission spectra [77] . Recent full-potential 
RPA calculations based on the maximally localized Wannier functions yielded U = 3.0 
eV [75], which is still too small in order to explain the photoemission spectra in terms 
of conventional DMFT calculations for the t2g band. This maybe a serious problem 
indicating that something is missing in the current interpretation of the photoemission 
data. Some of the missing ingredients may be the spacial correlations, the explicit 
contribution of the oxygen states, or the frequency- dependence of the effective Coulomb 
interaction in RPA [39j. On the other hand, the obtained value of the exchange 
interaction J' = 0.64 eV is very close to J' = 0.68 eV, which is typically used in 
the analysis of the photoemission spectra [76] . 
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In this small section we would like to consider two examples of deformation of the 
ideal perovskite structure, which are related with violation of the inversion symmetry 
around transition-metal sites. One is the oxygen vacancy, and the other one is the (001) 
surface of SrTiOs. Particularly, we will argue that such an inversion-symmetry breaking 
may lead to a number of new effects, and qualitatively change the character of transfer 
integrals between Wannier orbitals. 

6.2.1. Oxygen Vacancy in SrTiO^. In cubic perovskites, such as SrTiOs, the oxygen 
vacancy creates a dimer of Ti atoms located in its first coordination sphere. It also 
donates two electrons into the t2g bando 

In order to study the effect of the oxygen vacancy on the electronic structure of 
SrTiOa we have used the 3x3x3 supercell, in which one of the oxygen atoms has been 
replaced by the empty sphere. Such a composition corresponds to the chemical formula 
SrTi02.963- No lattice relaxation has been considered at this stage. According to LDA 
calculations, the electronic structure of such a supercell near the Fermi level is formed 
by 83 bands, which are well isolated from the rest of the spectrum (Figure [8]). Among 
them, 3x27=81 bands are the regular t2g bands, whereas two additional bands are 
formed predominantly by z^ orbitals of two Ti atoms located near the oxygen vacancy. 
The t2g and z^ bands are strongly mixed. 
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Figure 8. Left panel: LDA band structure of supercell Sr27Ti2708o corresponding to 
the oxygen-deficient SrTi02.963- The (dark) brown curves show 83 bands located near 
the Fermi level and formed by 3x27=81 t2g Wannier orbitals of all 27 Ti atoms as well 
as the z^ orbitals of two Ti atoms located near the oxygen vacancy (denoted as '1' 
and '2' in the right panel). The position of these z^ bands in the F-point are marked 
by the crosses. Right panel explains the atomic positions around the oxygen vacancy, 
namely, in the yz plane, which contains the vacancy, and in the xy plane, which is 
located just above the vacancy. The Ti atoms are indicated by the big dark (brown) 
spheres, the oxygen atoms are indicated by the small light grey (green) spheres, and 
the oxygen vacancy is indicated the symbol x . 



Therefore, it is clear that the minimal model near the Fermi level should be 

^■^ Under certain conditions, such a situation may lead to the formation of the spin-singlet bipolaronic 
state [71]. 
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constructed in the basis of four Wannier orbitals (nominally, z^, xy, yz, and zx) of 
two Ti atoms located near the vacancy, and three Wannier orbitals (nominally, xy, 
yz, and zx) of all remaining Ti atoms located in the next coordination spheres. The 
atomic wavefunctions of these types can be used as the trial functions in the downfolding 
method. The behavior of transfer integrals and the crystal-field splitting obtained after 
the downfolding is explained in Table [2l There is a number of interesting effects related 



Table 2. Crystal-field splitting (1-1) and transfer integrals in the bonds 1-2, 1-3, 1-4, 
and 1-5 around the oxygen vacancy in SrTiOs. All energies are measured in meV. The 
atomic positions are explained in Figure ID 
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with the presence of the oxygen vacancy. 

(i) The oxygen vacancy breaks the cubic symmetry and splits the t2g levels of two Ti 
atoms located next to it. The splitting is about 270 meV. However, already in 
the next coordination sphere, the t2g-level splitting is greatly reducedo and the 
situation becomes close to the perfect cubic environment. On the other hand, the 
position of the impurity z'^ level is lowered due to the missing Ti-0 bond. As a 
result, the z"^ levels become close to the t2g ones^ For example, the atomic splitting 
between the z^ and xy levels is only 40 meVcI 

(ii) The behavior of transfer integrals across the vacancy (the bond 1-2) is 
fundamentally different from the conventional case, when they are mediated by 
the oxygen 2p states (for example, in the bond 1-4): the transfer integrals between 
all three t2g orbitals are negligibly small, while the main interaction occurs between 
z'^ orbitals. 

^^ For example, for the titanium atoms 3 and 4 depicted in Figure [SJ the t2g-leve\ splitting is only 37 

meV and 40 meV, respectively. 

^® Note that the impurity z'^ level is an atibonding combination of the atomic oxigen 2p and titanium 

3z^ orbitals. Therefore, the lack of one of the Ti-0 bond formed by the Ti atom near the vacancy will 

shift the z^ level to the low-energy region. 

^^ For comparison, the t2g-eg splitting in the perfect perovskites is about 3 eV (Figures [1] and |3]). 
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(iii) The lack of the inversion symmetry leads to the mixing of the atomic 3d and Ap 
orbitals at the same Ti site. For example, the Wannier function, which is nominally 
denoted as z"^, besides the conventional Sd^i atomic orbitals will have some weight 
of the 4pz orbitals. Since the 4p orbitals are rather extended in the real space, 
such a mixing may change the form of the transfer integrals and even lead to the 
appearence of new interactions. The most striking example is the large transfer 
integral occuring between neighboring z^ and zx Wanner orbitals in the bond 1-3 
near the vacancy (Figure [9]). Such an interaction would vanishe in the perfect cubic 
environment. 
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Figure 9. Cartoon picture explaining the appearance of the "forbidden hoppings" 
near the points of the inversion symmetry breaking. The notations of the atomic sites 
are explained in Figure [5] Due to the local inversion-symmetry breaking (or the parity 
violation) , the atomic orbital d^^ at the site 1 is allowed to mix with the atomic orbital 
Pz, belonging to the same site. The combination of these orbitals will form a true 
Wannier function, centeerd at the site 1. Then, the electron can transfer from the Pz 
orbital of the site 1 to the dzx orbital of the site 3 (and vice versa). 

The distortion of the perfect cubic environment by the oxygen vacancy will affect 
not only the one-electron part of the model Hamiltonian (j3]), but also the Coulomb 
interactions (Table |3]). For example, around the vacancy, the Coulomb interactions 

Table 3. Intraorbital Coulomb interactions associated with different Ti sites around 
the oxygen vacancy in SrTiOs. All energies are measured in eV. The atomic positions 
are explained in Figure |S1 
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associated with different t2g orbitals are clearly different. This effect is captured by the 
RPA screening. The cubic symmetry of Coulomb interactions is practically restored is 
the fourth coordination sphere (site 5 in Figure [H]). 
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6.2.2. Surface states in SrTiO^. Another example of the inversion symmetry breaking 
is the Ti02 terminated surface of SrTiOa. The surface states has been studied in the 
slab geometry. Each slab contained nine Ti02 layers, which were separated by the SrO 
layers. Hence, the chemical formula of the slab was (Ti02)9(Sr0)8O According to the 
adopted notations, the first Ti02 layer corresponds to the surface, while the fifth Ti02 
layer corresponds to the bulk of SrTiOs (Figure [TOj) . No lattice relaxation has been 
considered at this stage. 
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Figure 10. Left panel: LDA band structure of the slab (Ti02)9(SrO)8. The (dark) 
brown curves are the t2g bands. The "surface" z^ bands are depicted by the symbols 
X. Right panel: The atomic positions near the (Ti02)-terminated surface of SrTiOa. 
The Ti atoms are indicated by the big dark (brown) spheres, and the oxygen atoms 
are indicated by the small light grey (green) spheres. 



The electronic structure of (Ti02)9(SrO)8 near the Fermi level consists of the 
3x9=27 t^g bands and two z^ bands, which are mainly formed by the surface Ti02 layers. 
Therefore, the minimal model can be constructed in the Wannier basis of (nominally) 
z^, xy, yz, and zx orbitals centered at the surface Ti sites and the xy, yz, and zx 
orbitals representing the remaining ( "bulk" ) sites. Thus, there is a direct analogy with 
the case of the oxygen vacancy in SrTiOa. In both cases, the Ti atoms located next to 
the "defect" (the surface, in the present case) acquire an additional z^ orbital, whereas 
other Ti sites are described in the standard t2g basis. 

The one-electron part of the model Hamiltonian is explained in Table HI At the 
surface, there is a huge crystal-field splitting, which even exceeds the crystal-field 
splitting near the single oxygen vacancy. Due to the inversion-symmetry breaking, 
there is an appreciable "forbidden" hopping between the z^ and zx Wannier orbitals 
operating in the surface bond 1-1'. The transfer integrals operating between t2g orbitals 
near the surface (the bonds 1-1' and 1-2) are also different from the ones in the bulk (the 
bond 5-5'). Thus, the effect of the surface on the electronic structure of the transition- 
metal perovskite oxides is not only in the narrowing of the yz and zx bands, caused 

^* The considered geometry can been obtained from the bulk SrTiOa by cutting the slab (Ti02)9(SrO)8 
and replacing the next (Ti02)2(SrO)3 layers by empty spheres. The considered region of empty spheres 
was sufficient to make the interaction between different slabs negligibly small. 



Combining DFT and Many-Body Methods to Understand Correlated Materials 



34 



Table 4. Crystal- field splitting (1-1) and transfer integrals at the surface and in the 
bulk of SrTiOa. All energies are measured in meV. The atomic positions are explained 
in Figure [TOl 
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by the reduced number of bonds available for the hoppings [63] O Even more serious 

consequences can be caused by the crystal-field splitting and the "forbidden hoppings" |fj 

The Coulomb interactions at the surface of SrTiOs are also considerably distorted 

in comparison with the bulk (Table |5]). However, the bulk-like behavior is practically 



Table 5. Intraorbital Coulomb interactions associated with different Ti-sites in the 
slab (Ti02)9(SrO)8. All energies are measured in eV. The atomic positions are 
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restored already in the second Ti02 layer. Similar to the single oxygen vacancy, the 
surface breaks the cubic symmetry of the Coulomb interactions. Moreover, the Coulomb 
interactions between t2g orbitals are somewhat smaller at the surface of SrTiOa then in 
the bulk. Like in the case of the single oxygen vacancy, this dependence of the effective 
Coulomb interactions on the local environment of the transition-metal sites is captured 
by the RPA screening. 

^^ Note that the orbital zx is perpendicular to the surface. Therefore, it can be involved in the hoppings 
in the directions ±a; and —z (in the geometry shown in Figure [T0|). Similar situation holds for the yz 
orbitals. On the contrary, the xy orbital is involved in the hoppings in all four directions ±a; and iy. 
^° Note that in addition to the hopping, the z^ and zx orbitals are coupled at the same transition-metal 
site by the spin-orbit interactions. If the surface were magnetic, this type of coupling would lead to the 
Dzyaloshinsky-Moriya interactions between the spins [801 [H]. Thus, the "forbidden hoppings" provide 
a microscopic basis for the appearence of these inetractions. 
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6.3. Distorted Perovskite Oxides 

The transition-metal perovskite oxides ABO3 (where A= Y or La, and B= Ti or V) 
are regarded as some of the key materials for understanding the strong coupling among 
spin, orbital, and lattice degrees of freedom in correlated electron systems [H |82] . 

All compounds have distorted perovskite structure. The distortion can be either 
orthorhombic (LaTiOa, YTiOa, and YVO3 below 77 K) or monoclinic (LaVOa and YVO3 
above 77 K). The space group of the orthorhombic phase is D^f^ (in Schonflies notations 
or Pbnm in the Hermann- Maguin notations. No. 62 in the International Tables). The 
space group of the monoclinic phase is Cf^^ (P2i/a, No. 14 in the International Tables) o 
The magnitude of the distortion is controlled by the size of the A^~^ ions. Typically, the 
smaller is the size of the A^~^ ions, the larger is the distortion. For example, the YBO^ 
oxides are always more distorted than the LaBO^ oxideso 

According to electronic structure calculations in the local-density approximation, all 
compounds have a common transition-metal t2g band, located near the Fermi level and 
well isolated from the rest of the spectrum (Figured]). The number of electrons donated 
by each transition- metal site into the t2g band is one and two for ATiO^ and AVO3, 
respectively. These electrons are subjected to strong Coulomb correlations, beyond the 
local-density approximation, and the systems are classified as Mott insulators [1]. 

Furthermore, the Coulomb correlations interplay with the lattice distortions, 
leading to a peculiar phase diagram, where each compound has a distinct magnetic 
structure (Figure [TTj) . For example, YTiOa is a ferromagnet [85l IMl [87l EH]. LaTiOa 
is a three-dimensional (G-type) antiferromagnet [81] 189] . At the low temperature, 
YVO3 has G-type antiferromagnetic structure, which transforms to a chain-like (C-type) 
antiferromagnetic structure at around 77 K (83] [901 EH |92]. The magnetic transition 
coincides with the structural orthorhombic-to-monoclinic transition. On the other hand, 
LaVOs is the C-type antiferromagnet in the whole temperature range below the magnetic 
transition point (93] [9l]. The understanding of these magnetic properties will be the 
main goal of this section. 

It is worth noticing that on the theoretical side the magnetic behavior of these "^29" 
compounds has been and continues to be the subject of numerous controversies. So far, 
the theoretical interpretation of the unusual magnetic properties has been based on two 
different standpoints, which practically exclude each other. One is the picture of the 
orbital liquid, which implies that the effect of the crystal distortion on the electronic 
structure is small and the (quasi-) degeneracy of the atomic t2g levels is preserved even 
in the real crystal environment (95] |96] . Another one is the theory of crystal field, which 
lifts the orbital degeneracy and stabilizes certain orbital structure being compatible with 
the observed magnetic ground state [251 EH 198] . 

•^^ There are certain indications that the actual symmetry can be even lower than P2i/a |83) . 
32 For example, the Ti-O-Ti angle is 157-158° in LaTiOg [HI], and only 140-144° in YTiOg [SS]. 
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Figure 11. Crystal and magnetic structure of distorted perovskite oxides. Tire La 
and Y atoms are indicated by the big (blue) dark spheres, the Ti and V atoms are 
indicated by the medium (red) dark grey spheres, and the oxygen atoms are indicated 
by the small (green) light grey spheres. The vectors a, b, and c show the directions 
of orthorhombic translations. The directions of the magnetic moments are shown by 
arrows. 



6.3.1. Results of Downfolding for the One- Electron Part and Wannier Functions. 
An example of the electronic structure in the region of the t2g bands for YTiOa is 
shown in Figure [121 It reveals an excellent agreement between results of the original 



YTiO, 




r xs Y rzu rt 



Figure 12. LDA bands for YTiOs as obtained in the original electronic structure 
calculations using the LMTO method and after the tight-binding (TB) parametrization 
using the downfolding method |22j . Notations of the high-symmetry points of the 
Brillouin zone are taken from 1651. 



LMTO calculations and their tight-binding parametrization obtained in the downfolding 
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method. Since cubic t2g and Cg orbitals are mixed by the crystal distortion, the trial 
functions {xr} of the downfolding method cannot be longer chosen from geometrical 
considerations. Generally, such a choice is not unique. Nevertheless, for these purposes 
one can use the eigenvectors obtained from the diagonalization of the density matrix. 
As it was already pointed out in Section 1331 such a procedure should guarantee a good 
degree of localization of the Wannier functions as well as of the parameters of the model 
Hamiltonian in the real space. 

The spacial extension of the Wannier functions constructed for YTiOs and LaTiOa 
is illustrated in Figure [131 In LaTiOs, about 80-85% of the total weight of the 



1.00 
0.95 

s 

.2? 0.90 



_ >^ 

S 0.85 I- D- 



0.80 
0.75 




Ti 



12 3 4 5 

distance from the central site (A) 




12 3 4 5 6 

distance from the central site (A) 



Figure 13. Spacial extension of Wannier functions for YTiOa (left) and LaTiOs 
(right). The figure shows the total weight of the Wannier function accumulated around 
the central Ti site after adding every new sphere of neighboring atomic sites. The 
symbols '1', '2', and '3' denote the Wannier functions obtained after projection onto 
the eigenvectors of the density matrix, where '1', '2', and '3' correspond to the largest, 
second largest, and third largest eigenvalues, respectively. 



Wannier function is accumulated at the central Ti site, 5-9 % belong to six oxygen 
sites surrounding the central Ti atom, and about 10 % is distributed over La, Ti, and O 
sites located in next coordination spheres. In YTiOs, the same distribution parameters 
are 82-87%, 6-10%, and 5%, correspondingly for the central Ti site, its neighboring 
oxygen sites, and Y, Ti, and O sites located in the next coordination spheres. Another 
measure of localization is the expectation value of square of the position operator: 
{r^)a={W^\{r-Rf\W^) fT6], which yields {r^)a= 2.68, 2.36, and 2.37 A^ for LaTiOa, 
and {r^)a= 2.28, 1.90, and 2.05 A^ for YTiOg. Thus, the Wannier functions for LaTiOg 
and SrVOs are less localized in comparison with the more distorted YTiOs. However, 
this is to be expected. One reason is the lattice distortion, which increases in the 
direction SrVOa — )■ LaTiOa — ?■ YTiOs. Another reason is the proximity of transition- 
metal t2g and oxygen 2p bandsO The degree of the localization in SrVOs and LaTiOa 
is very similar, despite the distortion which takes place in LaTiOs. One reason may be 

•^^ For example, the distance between the transition-metal t2g and oxygen 2p bands is 3.2 eV in YTiOs, 
2.7 eV in LaTiOs, and only 0.3 eV in SrVOa. Therefore in SrVOa, the mixing between the t2g and 
2p states is stronger, and the Wannier functions have a larger weight at the oxygen sites. To smaller 
extent the same is true for LaTiOa. 
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the proximity of the lanthanum 5d states to the Fermi level, which leads to stronger 
hybridization with the transition- metal t2g states |5| 199]. 

6.3.2. Transfer Integrals and Crystal-Field Splitting. The behavior of parameters of 
the one-electron part of the model Hamiltonian is explained in Figure [TH Because of 
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Figure 14. (a) : Crystal structure of the distorted perovskite oxides with the notation 
of the positions of four transition-metal sites in the unit cell. In orthorhombic 
compounds all sites are equivalent and can be transformed to each other by the 
symmetry operations of the 1)2^ gi'oup. In monoclinic compounds, there are two 
inequivalent pairs of the transition-metal sites: (1,2) and (3,4). (b): Distance 
dependence of averaged transfer integrals. In the orthorhombic (o) structure, all 
sublattices are equivalent and shown by a single symbol. In the monoclinic (m) 
structure, the transfer integrals around two inequivalent transition-metal sites are 
shown by closed and open symbols, (c): i2g-level splitting. The notations 'site 1' and 
'site 3' stand for two nonequivalent transition- metal sites in the monoclinic structure 
0. 



the lattice distortion, there is an appreciable crystal-field splitting, which is larger for 
the more distorted YTiOa and YVO3, and substantially smaller for the least distorted 
LaTiOs- Furthermore, there is a clear correlation between the number of t2g electrons 
and the form of the crystal-field splitting. As the rule, the crystal-field splitting tends 
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to quench the orbital degrees of freedom. For example, in YTiOs, it splits off one t2g 
level to the low-energy part of the spectrum, that is just enough to accommodate one t2g 
electron. The gap separating the lowest t2g level from the middle one is about 109 meV. 
On the contrary, the crystal-field splitting in YVO3 lowers the energies of simultaneously 
two t2g levels, that is again consistent with the number of t2g electrons per one vanadium 
site. The distance between middle and highest t2g levels in the orthorhombic phase of 
YVO3 is about 111 meV. 

The monoclinic distortion in YVO3 and LaVOa creates two inequivalent types of 
vanadium atoms, which lie in different ab-planes and are denoted as (1,2) and (3,4) 
in Figure [TH Typically, the less distorted planes (1,2) are alternated with the more 
distorted planes (3,4). The magnitude of the crystal-field splitting depends on the 
compound. For example, in YVO3, the energy splitting between middle and highest t2g 
levels in the planes (1,2) and (3,4) is rather similar: 101 meV and 128 meV, respectively. 
However, in LaVOa, the plane (1,2) appears to be much less distorted in comparison 
with the plane (3,4), resulting in different crystal-field splittings: 78 meV and 152 meV, 
respectively. Not only the energies, but also the directions of the crystal-field splitting 
are different in the orthorhombic and monoclinic phases, that immediately follows from 
the symmetry considerationso In Section I6.3.4[ we will argue that such a difference is 
directly related with the type of the orbital ordering, which is reflected in the magnetic 
properties of the distorted perovskite oxides. It is important to note that nonsphericity 
of the Madelung potential flT7|) is crucial for reproducing the correct magnetic ground 
state of YVO3 and LaV03, contrary to the conventional atomic-spheres approximation 

Because of the complexity of transfer integrals in the distorted perovskite structure, 
it is practically impossible to discuss the behavior of individual matrix elements of 
II^RR'll' Nevertheless, some useful information can be obtained from the analysis of 
averaged parameters 
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where d is the distance between transition-metal sites R and R'. All transfer integrals 
are well localized and practically restricted by the nearest neighbors, located at around 
4A (Figure [HJo). As expected, the transfer integrals between the nearest neighbors are 
generally larger for the less distorted LaTi03 and LaV03, and smaller for the more 
distorted YTi03 and YVO3. 

6.3.3. Effective Interactions. Matrix elements of the effective Coulomb interaction in 
the t2g band can be computed by applying the combined constrained DFT plus RPA 
approach, which has been explained in details in Section |H At each transition-metal 

•^■^ For example, in the orthorhombic phase, the sites 1 and 3 (Figure [Ti)) can be transformed to each 
other by the 180° rotations around the c-axis associated with the translation by c/2. In the monoclinic 
phase, such a symmetry operation is no longer available. 
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site, it yields the 3x3x3x3 [/-matrix of screened Coulomb interactions in the basis 
of Wannier orbit als. The intersite interactions are considerably smaller and can be 
neglected [22]. For the explanatory purposes, the [/-matrix can be fit in terms of 
two Kanamori parameters: the intraorbital Coulomb interaction U and the intraatomic 
exchange interaction J' [60] . The results of such fitting are shown in Table El There 



Table 6. Results of fitting of the effective Coulomb interactions for the distorted 
perovskite oxides in terms of two Kanamori parameters: the intraorbital Coulomb 
interaction U and the exchange interaction J7 [5]. All energies are measured in eV. 
The symbols 'o' and 'm' stand for the orthorhombic and monoclinic phase, respectively. 
The positions of the transition-metal sites are explained in Figure 1141 Generally, the 
site '1' is located in less distorted environment while the site '3' is located in more 
distorted environment. 
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is certain dependence of the parameter U on the local environment in solid, which is 
captured by RPA calculations [22]. For example, the value of W appears to be larger for 
the more distorted YBO-^ compounds. There is also a clear correlation between the value 
of U and the magnitude of the local distortion around two inequivalent transition-metal 
sites in the monoclinic structure: the sites experiencing larger distortion (according to 
the magnitude of the crystal-field splitting in Figure fT^) have larger L{, and vice- versa. 
On the other hand, the exchange interaction J' is less sensitive to the local environment 
in solids. 

Some of these trends can be rationalized by considering some additional 
approximations for the screening of Coulomb interactions in the t2g band by the oxygen 



2p band ( [Appendix Bj ). Nevertheless, we would like to emphasize that it is quite possible 
that these arguments are applicable only for this particular class of compounds. In other 
systems, other channels of screening may prevail and the situation can be generally 
different. 

Thus, even the effective Coulomb interaction in the t2g band can be affected by the 
crystal distortion, and this effect is clearly captured by realistic calculations of the model 
parameters. Besides two Kanamori parameters, U and J', the crystal distortion affects 
the fine structure of the entire [/-matrix and may even change its symmetry. Generally, 
such a symmetry is no longer cubic. This effect is clearly seen in the atomic multiplet 
structure, corresponding to two interacting t2g electrons in the case of vanadates: 
because of the crystal distortions the levels '^Ti^, ^E^, and ^T2g are slightly split (Figure 
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[T^ . The effect is not particularly strong. However, it can probably contribute to 
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Figure 15. Atomic multiplet structure of two interacting t2g electrons in tlie "ground 
state" of vanadates. Tire positions of the transition-metal sites in the monoclinic 
structure are explained in Figure [Hj Generally, the sites '1' are located in less distorted 
environment, while the sites '3' are located in more distorted environment. 



some delicate magnetic properties of the distorted perovskite oxides, such as the orbital 
magnetization, magnetocrystalline anisotropy, Dzyaloshinsky-Moriya interactions, etc. 
Finally, it is important to note that the obtained parameters of on-site interactions 
are substantially smaller than the value of Coulomb integral W ~ 5 eV derived from the 
photoemission data [76], which is typically used in calculations based on the dynamical 
mean-field theory |99l I100[ llOlj . As it was already pointed out in Section I6.1.2[ 
the reason for such a difference is rather transparent and related with the additional 
screening coming from the oxygen 2p as well as the transition-metal Cg bands, which are 
not explicitly included into the low-energy model (|3]). The application of the dynamical 
mean-field theory with the parameters reported in Table [6] would apparently lead to a 
metallic behavior for all considered compounds |99], being in straight contrast with the 
experimental data. Apparently, this contradiction can be resolved by going beyond the 
single-site approximation and considering the spacial correlations [52] . 



6.3.4- Orbital Ordering and Magnetic Properties. It is convenient to start the 
discussion with the low-temperature orthorhombic phase of YVO3, which is one of 
the most distorted systems among the considered ones. The distribution of the t2g 
electron density around the vanadium sites (the so-called orbital ordering) obtained in 
the Hartree-Fock calculations is shown in Figure [161 The type of the orbital ordering 
can be identified as "C" , meaning that the orientation of the orbital clouds relative to 
each other is nearly orthogonal in the ab plane, and nearly parallel along the c axis, i.e. 
in the close analogy with the C-type antiferromagnetic spin ordering (Figure [TTj) . 
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Figure 16. Distribution of the electron density around vanadium sites as obtained 
in the Hartree-Fock calculations for the ferromagnetic (F) and A-, C-, and G-type 
antiferromagnetic alignment in the orthorhombic phase of YVO3 (T < 77 K) [5]. 
Different magnetic sublattices are shown by different colors. 



The orbital ordering is very robust, and practically does not depend on the type of 
the magnetic arrangement. This means that the orbital degrees of freedom are quenched 
by the crystal-field splitting and the superexchange processes, which may adjust the 
orbital ordering in order to minimize the total energy of the system ^0\, are considerably 
weaker. 

Generally, the orbital ordering of the C-type is compatible with the G-type 
antiferromagnetic spin ordering [76l I102[ 1103] . Therefore, it is not surprising that the 
experimentally observed G-type antiferromagnetic ground state of orthorhombic YVO3 
is reproduced already at the level of mean-field Hartree-Fock calculations (Table [7j). 
The conclusion is totally consistent with results of all-electron LDA+f/ calculations 
|103[ 1104] , and provides a transparent physical explanation for them. 

The correlation effects beyond the mean-field Hartree-Fock approximation play a 
very important role and additionally stabilize the G-type antiferromagnetic ground state 
(Table [7j). The values of correlation energy obtained in the second order perturbation 
theory are summarized in Figure [T71 One can clearly see that for the orthorhombic phase 
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Table 7. Total energies of ABO3 compounds in the ferromagnetic (F) and A-, C-, 
and G-type antiferromagnetic states measured in meV per one formula unit relative 
to the most stable state (marked by the symbol x) as obtained in the Hartree-Fock 
calculations and after taking into account the correlation energy in the second order 
perturbation theory (shown in the parenthesis). The correlation energy includes both 
on-site and intersite contributions. The symbols 'o' and 'm' stand for the orthorhombic 
and monoclinic phases, respectively. 
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Figure 17. Correlation energies in the distorted perovskite oxides as obtained in 
the second order perturbation theory starting from the Hartree-Fock ground state for 
the ferromagnetic (F) and A-, C-, and G-type antiferromagnetic alignment. Only on- 
site contributions are shown. The symbols 'o' and 'm' stand for the orthorhombic 
and monoclinic phases, respectively. The positions of the transition-metal sites in 
the monoclinic structure are explained in Figure [T31 Generally, the sites T' are 
located in less distorted environment, while the sites '3' are located in more distorted 
environment. 



of YVO3, the correlation energy is clearly the largest in the G-type antiferromagnetic 
state among all considered magnetic structures. 

The behavior of interatomic magnetic interactions has been discussed in [5]. In 
the case of orthorhombic YVO3, there is a good agreement with experimental data 
concerning both the form (nearly isotropic three-dimensional antiferromagnetic network) 
and the absolute values of these interactions |HH] • 

The transition to the higher-temperature monoclinic phase of YVO3 is associated 
with an abrupt change of the orbital ordering patters, which can be identified as the 
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"G-type" and corresponds to nearly orthogonal orientation of the orbital clouds both in 
the ab-plane and along the c-axis (Figure [T8|) . Therefore, it is reasonable to expect the 





Figure 18. Distribution of the electron density around vanadium sites as obtained 
in the Hartree-Fock calculations for the ferromagnetic (F) and A-, C-, and G-type 
antiferromagnetic alignment in the monoclinic phase of YVO3 (T > 77 K) ^Sj. Different 
magnetic sublattices are shown by different colors. The arrows point at the orbital 
clouds in the (1, 2) planes, which slightly change their shape depending on the magnetic 
ordering. On the contrary, the orbitals in the (3,4) planes are well quenched by the 
crystal distortion. 



emergence of the C-type antiferromagnetic structure, which is generally compatible with 
the orbital ordering of the G-type |102j . In this sense, there is a close analogy with the 
low-temperature orthorhombic phase, where the type of the magnetic ordering is strictly 
determined by the type of the orbital ordering. However, there is also a difference. In 
the monoclinic phase, there are two types of vanadium atoms, which are located in the 
planes (1,2) and (3,4), respectively. The orbital degrees of freedom in the (3,4) plane are 
almost rigidly quenched by the large crystal-field splitting, similar to the orthorhombic 
phase. However, in the (1,2) plane this quenching is not complete, leaving some room 
for the orbital fluctuations. This effect can be seen even visually in Figure [181 where the 
shape of the orbital clouds in the (1,2) plane slightly change depending on the magnetic 
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state in order to additionally minimize the energy of superexchange interactions [TU] . 

Nevertheless, the mean-field Hartree-Fock approach appears to be a good 
starting point also for the monoclinic phase of YVO3. It correctly reproduces 
the C-type antiferromagnetic ground state (Table [7]), in agreement with all-electron 
LDA+U calculations |104] . The correlation energy additionally stabilizes the C-type 
antiferromagnetic phase relative to other magnetic structures. It is also important that 
the orthorhombic-to-monoclinic transition in YVO3 is well refiected in the behavior of 
the correlation energy, which was the largest for the G-type antiferromagnetic state 
in the orthorhombic phase and becomes the largest for the C-type antiferromagnetic 
state in the monoclinic phase (Figure [T71) . In the monoclinic phase, there is certain 
dependence of the correlation energy on the local environment of the vanadium sites: 
generally, the site '1' with less distorted environment has larger correlation energy (and 
vise versa), though the effect is not particularly strong. 

Similar to the monoclinic phase of YVO3, the orbital ordering in LaVOs tends to 
stabilize the C-type antiferromagnetic ground state. The orbital degrees of freedom are 
quenched in the (3,4) plane and retain enough flexibility in the (1,2) plane. For example, 
each change of the magnetic state is accompanied by the substantial reconstruction of the 
orbital ordering patter in the (1,2) plane, which is shown by arrows in Figure [T9l Thus, 
the importance of the orbital fluctuations should rise in the direction orthorhombic 
YVO3 — )■ monoclinic YVO3 — )■ (monoclinic) LaV03. Apparently, the application of 
the mean-fleld Hartree-Fock theory to the latter compound is already rather critical. 
Nevertheless, it still provides a consistent explanation for the number of properties 
of LaV03. For example, the C-type antiferromagnetic ground state is successfully 
reproduced by the Hartree-Fock calculations and additionally stabilized by correlation 
effects treated in the second order perturbation theory (Table [7j). The absolute values 
of correlation energy are larger than in YVO3 (Figure [T7|) . but still substantially smaller 
than the crystal-fleld splitting. Apparently, this is one of the reasons why the Hartree- 
Fock theory is still applicable. 

The situation with YTi03 is similar to the orthorhombic YVO3. The orbital degrees 
of freedom are frozen in some particular conflguration by the crystal-fleld splitting. This 
orbital conflguration practically does not depend on the magnetic state (Figure [20l). The 
type of the orbital ordering is compatible with the ferromagnetic ground state, which 
can be reproduced at the Hartree-Fock level. The same conclusion has been drawn in all- 
electron LDA+t/ calculations |103l 1105] . The correlation effects tend to destabilize the 
ferromagnetic ground state (Table [7] and Figure [T7|) . Although the ferromagnetic state 
has the lowest energy even after taking into account the correlation effects, the energy 
separation from the next A-type antiferromagnetic state is very fragileo However, this 
seems to be consistent with the relatively low Curie temperature {Tc ~ 30 K |87] ) 
observed in YTi03 [5J. 

•^^ More precisely, the on-site correlations tends to destabilize the ferromagnetic ground state, as it is 
clearly seen from Figure 1171 Small intersite correlations partially compensate this trend, again in the 
favor of the ferromagnetic alignment ;61] . 
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Figure 19. Distribution of the electron density around vanadium sites as obtained 
in the Hartree-Fock calculations for the ferromagnetic (F) and A-, C-, and G-type 
antiferromagnetic alignment in LaVOa [S]. Different magnetic sublattices are shown 
by different colors. The arrows point at the orbital clouds in the (f , 2) planes, which 
change their shape depending on the magnetic ordering. On the contrary, the orbitals 
in the (3, 4) planes are well quenched by the crystal distortion. 



LaTiOa is clearly an exception. The orbital degrees of freedom are rather flexible 
and not completely quenched by the crystal distortion (Figure [2T|) . The mean-field 
theory yields an incorrect magnetic ground state (antiferromagnetic A-type, instead 
of the G-type, Table [7]). The correlation interactions treated as a perturbation 
to the Hartree-Fock solution do not change the overall picture, and the G-type 
antiferromagnetic state remain unstable relative to the A-stateO 

Thus, the origin of the G-type antiferromagnetic ground state in LaTiOs remains 
a challenging problem for the future theories. Apparently, one of the basic assumptions 
about nondegeneracy of the Hartree-Fock ground state breaks down in the case of 
LaTiOs, and the true magnetic ground state cannot be continuously approached through 
the series of corrections applied to the Hartree-Fock ground state. One intuitive reason 

^^ Effects of the higher-order perturbation theory expansion have been considered in |6f| . in the 
framework of the T-matrix theory. Despite some quantitative differences from the second-order 
perturbation theory, the main conclusion was the same and the G-type antiferromagnetic state was 
always higher in energy than the A-type antiferromagnetic state. 
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Figure 20. Distribution of tlic electron density around titanium sites as obtained 
in the Hartree-Foek calculations for the ferromagnetic (F) and A-, C-, and G-type 
antiferromagnetic alignment in YTiOa [5,. Different magnetic sublattices are shown 
by different colors. 



could be seen already in Figure [T71 the correlation energies, derived as perturbations 
to the Hartree-Fock solutions, are exceptionally large in the case of LaTiOs, so that 
some of them become comparable with the crystal- field splitting (Figure [T^ . This 
means that there is certain inconsistency in our theoretical treatment: the perturbation 
becomes comparable with the basic energy splitting, that was used in order to justify 
this perturbation-theory expansion. 

In summary, the crystal-field theory can indeed explain many properties of the 
distorted perovskite oxides, but apparently not all. Particularly, the decrease of 
the crystal distortion in the case of lanthanides leaves some room for the orbital 
fluctuations. LaTiOs is definitely the most difficult case, where the crystal-field theory, 
if supplemented with realistic values of the model parameters, fails to reproduce the 
correct magnetic ground state. 

The main conclusions of this section are based on the second-order perturbation 
theory for the correlation interactions, which may be questionable because the values 
of the Coulomb repulsion U (Table [6]) and not particularly small and it is reasonable 
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Figure 21. Distribution of the electron density around titanium sites as obtained 
in the Hartree-Fock calculations for the ferromagnetic (F) and A-, C-, and G-type 
antiferromagnetic alignment in orthorhombic phase of LaTiOa [5] . Different magnetic 
sublattices are shown by different colors. 



to expect that higher-order effects may modify some of these conclusions. The role 
of higher-orders of the perturbation theory expansion for the correlation energy was 
investigated in [61], in the frameworks of the T- matrix theory. Generally, the higher- 
order effects tend to reduce the correlation energy. However, the amount of this 
reduction strongly depends on the magnetic state. For example, the correlation energy 
in the ferromagnetic state is overestimated in the second-order perturbation theory by 
less than 15% in comparison with the T-matrix method. For the antiferromagnetic 
configurations, the effect is more pronounced and the correlation energy can be 
overestimated by about 50%. Nevertheless, the main trends in the behavior of the 
correlation energy are well captured already by the second-order perturbation theory. 
For example, if the absolute value of the correlation energy was the largest for certain 
magnetic configuration in the second-order perturbation theory, the same tendency 
is clearly seen in the T-matrix theory, etc. In this sense we believe that the main 
conclusions of the present Section are valid. Nevertheless, the correlation interactions 
in the distorted perovskite oxides is certainly one of the most important and interesting 
problems, which deserves a thorough investigation. The analysis based on the second- 
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order perturbation theory, considered in this Section, can be regarded as only the first 
step in this direction. 

6.3.5. Effects of Spin- Orbit Interaction. The spin-orbit interaction in distorted 
perovskite structures generally leads to a noncoUinear magnetic alignment, which obeys 
certain symmetry rules |106[ I107[ I108[ 1109] . The spin magnetic moments aligned along 
one of the orthorhombic axes are subjected to certain rotational forces originating 
from anisotropic and Dzyaloshinsky-Moriya interactions |80l [81], which lead to the 
reorientation (or canting) of these magnetic moments. Hence, in the equilibrium 
magnetic structure we will generally have all three projections of the magnetic moments 
onto the orthorhombic axes. Furthermore, the type of the magnetic ordering for these 
three projections will be generally different. Thus, each magnetic structure can be 
generally abbreviated as X-Y-Z, where X, Y, and Z is the type of the magnetic ordering 
(F, A, C, or G) formed by the projections of the spin magnetic moments onto the 
orthorhombic axes a, b, and c, respectively. The orbital magnetic structure has the 
same sjTiimetry, although it may have a different origin of the canting, which arises 
mainly from the minimization of the single-ion anisotropy energy at each transition- 
metal site. Generally, the spin and orbital magnetic moments are not coUinear to each 

other [TnmiTTn] . 

For example, the true magnetic ground state of YTiOs is G-A-F (Figure |22|) . 
where in addition to the main ferromagnetic components along the c axis, there will 
be two antiferromagnetic components parallel to the a and b axes. In this case, one 
can expect rather interesting consequence of the correlation effects on the noncoUinear 
magnetic alignment [5]. As it was already pointed out in previous Section, without 
spin-orbit interaction, the correlation effects in YTiOa favor the antiferromagnetic 
interatomic coupling and systematically lower the energies of all antiferromagnetic 
structures relative to the ferromagnetic one. Then, it is reasonable to expect that 
after taking into consideration the spin-orbit interaction, the correlation effects will 
systematically increase the weight of the antiferromagnetic components in the ground 
state configuration and result in the additional spin canting away from the collinear 
ferromagnetic state. Such an effect is clearly seen in Figure [221 after taking into account 
the correlation effects in the frameworks of variational superexchange theory, about 8% 
of the spin magnetization-density is transferred from the ferromagnetic part parallel to 
the c axis to the antiferromagnetic part lying in the ab plane. The distribution of the 
orbital magnetization-density is less sensitive to the correlation effects. We also note 
that the correlation effects readily explain the experimental values reported for the F 
and G components of the magnetic moments j88j . 

Another interesting phenomenon, which is directly related with the spin-orbit 
interaction, is the temperature-induced magnetization reversal behavior observed in 
YVO3 (Figure [23]) [90j . Briefly, upon cooling, the magnetization parallel to the a axis 
(in the notations of the Pbnm group) changes the sign two times: first continuously, 
in the monoclinic regions, and then discontinuously, at the point of monoclinic-to- 
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Figure 22. Spin (a and c) and orbital (b and d) magnetic structures of YTiOa in the 
Hartrce-Fock approximation (a and b) and after taking into account the correlation 
effects in the frameworks of variational superexchange theory (c and d) . The directions 
of the magnetic moments at different titanium sites are shown by arrows. The size of 
each arrow is proportional to the size of the magnetic moment. For the clarity, the 
arrows for the orbital magnetic moments have been additionally scaled by the factor 
ten. Corresponding values of spin and orbital magnetic moments, measured in Bohr 
magnetons, are given in the parenthesis. 



orthorhombic transition. Although the discontinuity of the magnetization is certainly 
expected for the first order nionoclinic-to-orthorhombic transition, its behavior in the 
monoclinic phase is somewhat surprising. Nevertheless, it is important to note that a 
similar effect is know to occur in some ferrimagnetic materials, like C02VO4, consisting 
of two or more types of antiferromagnetically ordered magnetic ions |lllj . Although 
there are several additional complications in comparison with C02VO4, this may be a 
clue for understanding the unusual behavior of YVO3 [5]. 

(i) As it was already pointed out in Section I6.3.4| although all vanadium atoms 
are chemically equivalent, in the monoclinic phase they are located in different 
crystallographic environments and one can clearly distinguish two types of 
vanadium atoms lying in the planes (1,2) and (3,4) (Figure [IE]). The magnetic 
interactions in these two planes are also different: while J34 is strongly 
antiferromagnetic, similar to the orthorhombic phase, J12 is considerably weaker^ 
Thus, there is a clear analogy with two magnetic sublattices of Co and V existing 
in C02VO4. However, in the monoclinic phase of YVO3, this difference is entirely 

•^^ Note that the existence of two magnetic sublattices also leads to a splitting of the magnon spectrum 
into acoustic and optical branches 1104] , which was clearly seen in the experiment [92j . 
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Figure 23. Temperature dependence of the magnetization in YVO3 (from [90]). The 
magnetic field of 1 kOe is applied along the a, b, and c axes, respectively. In the 
monoclinic phase, realized between 116 K and 77 K, the magnetization after first 
increasing start decreasing and crosses zero at around 95 K to a large negative value. 
With further cooling, it jumps to a large positive value at the point of monoclinic-to- 
orthorhombic transition. The orientation of the a, b, and c axes corresponds to the 
Pbnm group, both in orthorhombic and monoclinic phases. 

Table 8. Interatomic magnetic interactions (in meV) in the orthorhombic and 
monoclinic phases of YVO3, as obtained in the Hartree-Fock calculations without 
spin-orbit interaction f5^. The orthorhombic phase is stabilized below 77 K, while the 
monoclinic phase is stabilized above 77 K. The atomic positions are explained in Figure 
1141 The magnetic interactions were computed using formula (1291) for infinitesimal 
magnetic rotations near the antiferromagnetic ground state of the G-type in the case 
of the orthorhombic phase and the C-type in the case of the monoclinic phase. 
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related with the local crystal distortions, which create two inequivalent types of 
vanadium sites. 

(ii) An additional complication comes from the fact that YVO3 is an antiferromagnet, 
and no net magnetic moment is expected for the neither C- nor G-type 
antiferromagnetic ordering, realized in the monoclinic and orthorhombic phases, 
respectively. Nevertheless, it is reasonable to expect a weak ferromagnetism 
arising from the spin-orbit interaction in the distorted perovskite structure. In 
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the orthorhombic phase, the weak ferromagnetic moment is indeed ahgned along 
the a axis (in the notations of the Pbmn group), that directly follows from the 
symmetry considerations (Figure 12^ . In the monoclinic phase, the situation is a 
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Figure 24. Spin (a and c) and orbital (b and d) magnetic structures realized in the 
orthorhombic (a and b) and monoclinic (c and d) phases of YVO3 in the Hartree-Fock 
approximation. The directions of the magnetic moments at different vanadium sites 
are shown by arrows. The size of each arrow is proportional to the size of the magnetic 
moment. For the clarity, the arrows for the orbital magnetic moments have been 
additionally scaled by the factor ten. Corresponding values of spin and orbital magnetic 
moments, which are measured in Bohr magnetons, are given in the parenthesis. The 
orientation of the a, b, and c axes corresponds to the Pbnm group in the orthorhombic 
phase and the P2i/a group in the monoclinic phases. The P2i/a notations can be 
transformed to the Pbnm notations by interchanging the axes a and b |91) . 



little bit complicated because there is a large cancelation for all three projections of 
the magnetic moments. However, since the planes (1,2) and (3,4) are inequivalent, 
this cancelation is not complete and the net magnetic moment can be expected 
along the b direction in the notations of the P2i/a group, which corresponds to 
the a direction in the notations of the Pbmn group. 

Thus, there is at least a qualitative consistency with the experimental data [90j • A 
quantitative theory of the magnetization reversal behavior in YVO3 is still missing, and 
would be an interesting step to do. 



Combining DFT and Many-Body Methods to Understand Correlated Materials 



53 



6.4- Hyperoxide KO2 

As a final example, we would like to consider a molecular analog of "strongly-correlated" 
system without any transition-metal or rare-earth elements. 

The simplest molecule, whose ground state is different from the conventional spin 
singlet is O2HJ If these molecules form a crystal (by either cooling or pressurizing), 
it may become magnetic. Such a situation is indeed realized in the solid oxygen 
|112[ I113J . Since O2 is a good oxidizer and can easily take an additional electron when 
it brought in contact with alkali elements, there is another way of making crystalline 
arrays of the oxygen molecules, in the form of ionic crystals. Such materials do exist. 
One typical example is KO2, which was used as the chemical oxygen generator in 
rebreathers, spacecrafts, and life support systems. The magnetic properties of KO2 and 
other alkali hyperoxides were intensively studied in 1970s |114l 11151 11161 111711118 . 1119] . 
However, many details of their behavior remain largely unknown, even todaycj KO2 
has six crystallographic modifications, apparently related with different orientations of 
the oxygen molecules (Figure 1251) . However, the details are known only for the body- 
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Figure 25. Brief summary of structural and magnetic properties of KO2 |115| I116J . 

The oxygen molecules are shown as black dimers and potassium atoms are indicated 
by white symbols, /icff is the effective magnetic moment and Tew is the Curie- 
Weiss temperature. Both are derived from the analysis of magnetic susceptibility 
data in different temperature regimes. The arrow shows the region of "ferromagnetic- 
antiferromagnetic crossover" , where the Curie- Weiss temperature changes the sign. 



38 



The O2 molecule has two unpaired electrons in the doubly degenerate Hg level, resulting in the 
spin-triplet ground state, as it is required by the first Hund rule. 

^^ Evidently, the early progress was severely hampered by many objective difficulties related with the 
preparation of rigidly held samples and separation of intrinsic properties of KO2 from inhomogeneity 
effects. 
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centered tetragonal (bet) phase stabilized between 231 K and 400 K. Below 7K, KO2 
develops a long-range antiferromagnetic order |115[ 1120] . However, its type is unknown. 
In the paramagnetic region, the Curie- Weiss temperature (Tew) changes the sign at 
around 230 K, indicating at some change of interatomic magnetic interactions. Formally, 
this change may be related with one of the structural phase transitions. However, no 
profound change of the magnetic susceptibility has been observed around this transition, 
indicating that there might be a more fundamental reason for the change of Tew- 
Finally, the effective magnetic moment in the paramagnetic region (fies) is about 2 fiB, 
meaning that the large orbital contribution (about 1 /i^) is not quenched by the crystal 
distortion and persists down to 12 K. This suggests that the actual crystal distortion 
related with the reorientation of the oxygen molecules near the bet phase is quite week, 
and the bet phase itself can be regarded as a good starting point for the analysis of the 
magnetic properties of KO2 in the wide paramagnetic region |121j . 

The oxygen molecule appears to be the building block of not only the crystal, 
but also of the electronic structure of KO2 in the local-density approximation. The 
strong hybridization within the molecule leads to the formation of the molecular levels. 
The interaction between the molecules is considerably weaker, so that the molecular 
orbitals form a group of narrow nonoverlapping bands (Figure [26]) . Thus, there is a 
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Figure 26. Crystal structure and LDA density of states of the hct phase of KO2. The 
shaded area shows the contributions of the oxygen 2p states. Other symbols show the 
positions of the main bands. Fermi level is at zero energy. 



clear analogy with the atomic limit in the physics of strongly-correlated systems [TOl |82] , 
except that now, the localized electrons (or holes) reside on the molecular orbitals, which 
are distributed between two atomic sites. The doubly-degenerate vr^ band located near 
the Fermi level is formed by antibonding molecular px and Py orbitals. In comparison 
with the solid oxygen, the potassium atom donates an extra electron into the vr^ band. 
Therefore, the band is | filled. Due to the peculiar | filling, not only spin but also orbital 
degrees of freedom appear to be active and contribute to spin and lattice dynamics of 
KO2. In this sense, KO2 can be regarded as a molecular analog of correlated electron 
systems, comprising of orbitally degenerate magnetic O^ ions |121] . 
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6.4- 1- Construction and Parameters of Model Hamiltonian. Obviously that the 
minimal model for KO2 should be constructed in the basis of TCg bands and take into 
account the Coulomb correlations beyond conventional LDA. The basic idea of the 
construction of such a model is to relate each lattice point to the single oxygen molecule 
and formulate the problem in the Wannier basis corresponding to antibonding molecular 
Px and Py orbitals. 

The transfer integrals operating between such molecular Wannier orbitals can be 
derived by using formal downfolding procedure described in Section 13.31 The behavior 
of these transfer integrals in the bet lattice is explained in Figure |271 Basically, all 




Figure 27. Transfer integrals (in meV) associated with different bonds in the bet 
phase of KO2 |121j . The antibonding molecular orbitals are shown in the right part of 
the figure. The order of orbitals is p^ and Py. 



interactions are restricted by four nearest neighbors. Other interactions are considerably 



smaller. The on-site part of the one-electron Hamiltonian 
relativistic spin-orbit interaction, 
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in the basis of |pxt)) |Pj/t)) |Pxi), and \py\) orbitals, where the arrows indicate the spin 
state. The parameter of spin-orbit interaction ^ is about 17 meV [121] . 

The screened Coulomb interactions in the Hg band are computed in two steps [22] . 
First, the interaction parameters between atomic 2p orbitals are derived by using the 
constrained DFT method. It yields the intraatomic Coulomb interaction u ~ 11.4 eV, 
the interatomic intramolecular Coulomb interaction f ^ 1.3 eV, and the intraatomic 
exchange urteracfou , « 2.3 eV0 After that, it is uecessary to consider the additional 

""^ For comparison, the parameters of bare interactions are u w 27.8 eV, u ~ 11.0 eV, and j sa 2.9 eV. 
Thus, similar to the cubic perovskites considered in Section [6.1.2l the interatomic Coulomb interactions 
V appears to be the most screened, while the intraatomic exchange interaction j is the least screened. 
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screening of these interactions in the Tig band by other bandsQ and derive the full 
matrix f/=||f/a/37<5|| of screened Coulomb interactions between antibonding molecular 
orbitals px and Py, which is fully specified by the parameters of intraorbital Coulomb 
interaction W ~ 3.6 eV and the exchange interaction J' p^ 0.6 eV. This part can be done 
in the random-phase approximation, by starting from interaction parameters obtained 
in the constrained DFT. It is also important to remember that the interorbital Coulomb 
interaction W is related with U and J' by the identity W = IA — 2J. 

6.4-2. Implications to the Properties of K02. Since the Coulomb interactions U and W 
are clearly the largest parameters in the problem, the Hubbard model ([3]) can be further 
converted into a spin-orbital superexchange model by starting from the limit of isolated 
oxygen molecules and treating all transfer integrals as a perturbation [10]. Due to the 



I filling, it is convenient to use the hole representation and ascribe to each molecular 
site a single hole spin-orbital 

I a) = ab^\px t) + ab^\Py t) + ab^lp^ i) + ab^\Py i), (34) 

where |ap-|-|ap = |&tP+l^tP ~ I^^P+I^-;P = 1- Then, the energy gain 7rr' = 
T(aR, q;r') caused by the virtual hoppings in the bond (RR') can be computed using the 
formula fl33l) . In this case, G = G'(q;r,q;r/) is the Slater determinant constructed from 
the hole orbitals at the molecular centers R and R', -Er'm and |R'M) are the eigenvalues 
and eigenstates of the excited two-hole configurations at the molecular center R', and 
^R' is a projector operator enforcing the Pauli principle and preventing any hoppings 
of holes into or' . 

The problem can be further simplified by eliminating the orbital degrees of freedom, 
described by the 6-vasriable in ( IMl) . and constructing an effective spin model separately 
for each temperature regime. This can be generally done by averaging (!33|) with some 
distribution function V{b'^,b^,bi,bi,T) [10]. The formulation is especially simple for 
two limiting cases: T — )■ 0, corresponding to an orbitally ordered state, and T — )■ oo, 
corresponding to the complete orbital disorder. 

In the low-temperature limit kBT<^^ (note that C,/kB ~ 200 K), V{b^,bf,bi,bi,T) 
is fully controlled by the relativistic spin-orbit interaction, which picks up a 
linear combination of two spin-orbitals, \p+t) = ~i\Pxt)~APyt))/V^ and |p-t) = 
(|Pa;t)+^|Pyt))/v^, minimizing the spin-orbit interaction energy. Since each hole-orbital 
a is confined in the two-dimensional subspace spanned by |p+t) and |p-t)) the energies 
(133|) can be further mapped onto the anisotropic Heisenberg model with pseudospin 1/2 
[121]: 

^Heis = -^ E {(^R^R' + ^R^R') -^RR' + ^R^R'4r'} , (35) 

RR' 

"^^ Namely, the tTu band, that is a bonding combination of the p^ and py orbitals, as well as the a* and 
cr* bands, constructed from the pz orbitals (Figure 
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where 





Ty = \ . „ 1 , and f,. 



are the Pauh matrices in the basis of |p+t) ^^id |p_i) orbitals. The parameters of 
superexchange interactions, J^^t and J^j^^, can be expressed through the pair-interaction 
energies Trr' = T{a-R,aii') by considering different combinations of apj, and aR/O 
They are summarized in Table M The second-neighbor couphng Jp{^R^_,_b stabihzes the 



Table 9. Parameters of Heisenberg model (in meV) for the bet phase of KO2 in 
different temperature regimes |121j . J^^, and J^j^/ correspond to the low-temperature 
limit, while Jrr' corresponds to the high-temperature limit. The vectors separating 
two oxygen molecules, b = R' — R, are denoted as bi=(0,a, 0), b2 = (f,|,|), 
h3—{a,a,0), and b4=(0, 0,c) (see Figure [26l for the notations). 



vector 


7II 


7-L 
'^RR' 


Jrtl' 


bi 


-0.43 


-1.15 


-0.47 


b2 


1.31 


-1.09 


-1.01 


b3 


-0.02 


-0.05 


-0.02 


b4 


-0.07 


-0.07 


-0.01 



easy-axis ferromagnetic state. Other interactions are antiferromagnetic and frustrated 
on the bet lattice. Then, the Curie temperature can be estimated using the renormalized 
spin-wave theory |122j . It yields Tc ~ 70 K, which does not seem to be fully consistent 
with the experimental data, because no long-range magnetic order has been observed 
down to the very low temperature |1141I115| . However, there is also a number of factors, 
which may further affect this theoretical picture |121] . One is the deformation of the 
relativistic spin-orbitals |p+t) and \p—l) caused by superexchange interactions. Another 
one is the reorientation of the oxygen molecules (which probably takes place but not 
precisely known from early experimental studies |115[ 1116] ). For example, one can 
suggest that above 12 K, the reorientation of the oxygen molecules is not particularly 
large and leads only to some quantitative redefinition of parameters of the anisotropic 
Heisenberg model [351 which systematically decreases Tc |121] . However, it seems that 
below 12 K the situation changes dramatically. Apparently, the transition to the new 
crystallographic phase is accompanied by the large reorientation of the oxygen molecules, 
which not only quenches the orbital magnetic moment, but also plays a decisive role in 
the formation of the long-range antiferromagnetic order |121] . The complete quantitative 
theory describing these reorientation effects is missing at the present stage but would 
certainly be an interesting step to develop in the future. 

42 Namely, by denoting as |t±) = i\p+t)±\p-i))/V2, |t±) = i\p+t)±i\P-i))/V2, |r+) - |p+t>, 
and \t^) = \p-l) the pseudospin states corresponding to the positive (-I-) and negative (— ) directions 
parallel to the x-, y-, and z-exes, the parameters of superexchange interactions can be found as 2J^j^, = 

Ti^L^ Tr'x) - Ti^L^ Tr'J = '^i^Ry^^R'y) " '^(^Ry' ^R'y)' and 2j|j^, = T{t+^,t^,,) - T{t+^,t+,J. 
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Nevertheless, the intrigue of KO2 is that there is another region of antiferromagnetic 
interactions, which is reahzed in the high-temperature bet phase, as it is clearly 
manifested in the behavior of inverse magnetic susceptibility |114l 1115] . The origin 
of this "high-temperature antiferromagnetism" may be directly related with the orbital 
disorder. Indeed, in the high-temperature limit ksT^^, the thermal fluctuations will 
eventually destroy the relativistic coupling between spin and orbital degrees of freedom. 
Therefore, it is reasonable to assume that all orbital configurations are realized with 
equal probabilities. In this case, the spin system becomes fully isotropic, and the 
parameters of the spin Hamiltonian, 2 Jrr/ = 7^^/ — T^r/, can be derived by averaging 
the energies ( l33l) for antiferromagnetic (ti) and ferromagnetic (tt) configurations of 
spins in each bond over all combinations of orbital variables p^ 1121] . The orbital 
disorder gives rise to antiferromagnetic interactions Jrr', which are also summarized 
in Table [9l The corresponding Curie- Weiss temperature can be estimated using 
renormalizes spin-wave theory. It yields Tew ~ ~80 K, which is comparable with 
the experimental data |114l 1115] . 

In summary, the magnetic properties of KO2 provide an example of spin-orbital 
superexchange physics realized in the molecular solid. The properties largely depend 
on the orbital state of the O^ ions. In the paramagnetic region, the character of 
intermolecular interactions gradually changes from mainly ferromagnetic, and driven 
by the relativistic spin-orbit interaction, to antiferromagnetic, and corresponding to 
the picture of an independent spin and orbital disorder. All these features seem to 
be consistent with the experimental behavior of KO2. Finally, it is important to note 
that the geometry of molecular orbitals can be used as an additional degree of freedom, 
which controls the properties of superexchange interactions. Definitely, it adds many 
new functionalities into the classical problem of superexchange, which deserve further 
exploration. 

7. Summary and Concluding Remarks 

The goal of this review article was to outline the main ideas and scopes of new developing 
direction for the strongly correlated systems which can be called as the "realistic 
modeling" . The primary purpose of this project is to make a bridge between methods of 
first-principles electronic structure calculations, based on the density-functional theory, 
and many-body models, describing properties of strongly correlated systems in terms 
of a limited number of the most relevant model parameters and including information 
about all remaining electronic structure implicitly, through the renormalization of these 
model parameters. 

The realistic modeling has all potentials to become a useful tool of electronic 
structure calculations for the strongly correlated systems. It is true that we still have 
to rely on a number of approximations, particularly in the process of calculation of 
screened Coulomb interactions, which are inevitable in any approach. However, it is very 
important that, apart from these approximations, the entire procedure is parameter-free. 
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Namely, we do not have to deal anymore with adjustable parameters and we do need to 
resolve numerous ambiguities with the choice of these parameters. Instead, the realistic 
modeling brings the state of the discussion to a qualitatively new level: how to improve 
the approximations used for the definition and calculation of these model parameters. 

One undeniable advantage of realistic modeling is that it allows us to combine the 
accuracy and predictable power of first-principles electronic structure calculations with 
flexibility and insights of the model analysis. This idea was illustrated on the series of 
examples, for which we were able to consider the whole chain of actions starting from 
conventional electronic structure calculations in the local-density approximation, which 
was followed by the construction of an appropriate low-energy model, motivated by these 
calculations, and finally - by the solution of this model and by the analysis of properties 
of strongly correlated systems in terms of these model categories and trends. The first 
applications are indeed very encouraging and we would like to hope that in future the 
ideas of realistic modeling will continue to develop in order to become a powerful tool for 
theoretical analysis, design, and control of properties of strongly correlated materials. 
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Appendix A. Properties of Coulomb Interactions in the Atomic Limit 

In the atomic limit, the 5x5x5x5 matrix u of interactions between d electrons is totally 
specified by three radial Slater's integral F°, F"^, and F^, which are related with the 
parameters u and j, obtained in the constrained DFT, by the identities 

M = F° (A.l) 

and 

J = ^iF" + F") ■ (A-2) 

In order to derive all three Slater's integrals by knowing only two parameters u and 
j, extracted from the constrained DFT, one can use the additional relation 

F^/F"^ ^ 0.63, (A.3) 

which holds approximately in the atomic limit. 

We would like to emphasize that the parameter j, defined as flA.2p . is the measure 



of the exchange splitting corresponding to the spherically averaged electron densities 
for the majority- and minority-spin states in LDA. In the literature one can find other 
definitions of the exchange integrals, which sometimes cause certain confusion. For 
example, Mizokawa and Fujimori [76] defined the exchange integral as 

jMF = \b + C, (A.4) 
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in terms of the two Racah parameters: B = {9F^ — 5F^)/441 and C = 5F^/63. Then, 
it is easy to verify the vahdity of the following relation between j and Jmf: 

JMF = ^J ~ 0.71J. (A.5) 

Thus, jmf is always smaller than j. 

Alternatively, one can define parameters of Coulomb interactions only for the t2g 
manifold. There are three types of interactions, which are sometimes called as the 
Kanamori parameters [60] : 

u, = J dr J dr'W^^{r)W,y{r)\r-rr'W^y{r')W.y{r'), (A.6) 

[ = j dv j dv'Wl{v)W,y{v)\v-v'\-'Wl{j')Wy,{r'), (A.7) 



u 



and 

it 



I dr I dr'Wl^{r)Wy.{r)\r -Tr'Wl^{r')Wy.{r'). {A.i 



In the atomic limit, they can be expressed in terms of radial Slater's integrals as 

„,.f<. + i.f= + _Lf.. (A.9) 

and 

jt = —F^ + —F\ (A.ll) 

•^ 49 441 ^ ' 

Other types of interactions between t2g electrons vanish. In fact, there are only two 

independent interactions because Ut, u[, and jt are related by the identity: 

ut = u[ + 2jt. (A.12) 

It is also straightforward to show that 

ut = u + -J. (A.13) 

The parameter jt can be expressed through j using the approximate relation between 
Slater's integrals (]A.3|) . which yields 

Jt ~ 0.77J. (A.14) 

Thus, generally we have the following inequality for the exchange integrals defined in 
three different ways: 

jMF<Jt<J- (A. 15) 

After taking into account the RPA screening, the parameters Ut, u[, and jt 
correspond to the parameters U, W, and J' considered in Sections 16.1.21 and 16.3.31 
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For the /-shell, the 7x7x7x7 matrix u can be reconstructed from F^, F"^, F'^, and 
F^, by using the following identities in the atomic limit [71] : 

u = F°, (A.16) 

7 = - (—F^ + —F^ + —fA (A 17) 

^ 3 Vl5 11 ^429 J' ^ ' 

pijF^ ^ 451/675^ and F^F^ ~ 1001/2025. 

Appendix B. Correlations Between 2p-3d Hybridization and Screening of 
Coulomb Interactions in the t2g band of transition-metal perovskite oxides 

In this appendix we derive some approximate expression for the static RPA screening of 
Coulomb interactions in the t2g band of transition-metal oxides by the oxygen 2p band. 

First, we assume that the band dispersion is considerably smaller than the charge- 
transfer energy Act? which is the energy difference between the centers of gravity of the 
oxygen 2p band and the transition-metal t2g and Cg bands (Figure[l]). Then, for the static 
{u = 0) screening caused by the oxygen 2p band, the denominator of the polarization 
function (|24l) can be replaced by Act and one can perform separate summation over 
the occupied and empty states. 

Then, we focus on the self- screening caused by the atomic 3d states, which 
contribute to the oxygen 2p band due to the hybridization effects and consider the 
matrix elements of (12^ in the basis of atomic 3d orbitals. These matrix elements can 
be expressed through the density matrices of the 3d states calculated separately in the 
occupied oxygen 2p band and in the empty part of the spectrum (correspondingly h° and 
n'^). Finally, we assume that n°f ~ N''^'^^5a(s {N° and N^ being the total number of 3d 
electrons in the oxygen 2p band and in the empty part of the spectrum, respectively) and 
consider the diagonal matrix elements of (12^ . which mainly contribute to the screening 
of the intraorbital Coulomb interaction U. In this case, one can to derive the following 
(approximate) expression for the diagonal matrix elements of the polarization: 

V 



Act 

Moreover, since the total number of 3d electrons is conserved, N°+N'^ is a constant, 
which depends only on the number of 3d electrons in the occupied part of the t2g band. 
Therefore, it is reasonable to expect the following rules for compounds having the same 
number of t2g electrons: 

(i) the larger is the weight of the transition-metal 3d states in the oxygen 2p band, A^", 
the stronger is the screening of Coulomb interactions in the t2g band; 

(ii) N° is controlled by the hybridization between oxygen 2p and transition-metal 3d 
states, which in turn depends on the crystal distortion (particularly, the buckling of 
the Ti-O-Ti and V-O-V bonds). Therefore, stronger distortion will tend to reduce 
the screening; 
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(iii) smaller Act in the case of SrVOs will additionally increase the screening of Coulomb 
interactions in the t2g band. 

All these trends are clearly seen in Table IBll which shows the correlation between 
N° and U in different transition-metal perovskite oxides. For example, more distorted 

Table Bl. Correlation between number of 3d electrons in the oxygen 2p band {N°) 
and the value of screened Coulomb interaction {U^ in eV) in the transition- metal <2g 
band. The symbols 'c', 'o', and 'm' stand for the cubic, orthorhombic, and monoclinic 
structure, respectively. The positions of the transition-metal sites are explained in 
Figure [131 Generally, the site '1' is located in less distorted environment while the site 
'3' is located in more distorted environment. 



compound 


phase 


site 


N° 


U 


YTiOa 







0.66 


3.45 


LaTiOa 







0.73 


3.20 


SrVOa 


c 




1.19 


2.53 


YVO3 







0.74 


3.27 


YVO3 


m 




0.76 


3.19 






3 


0.72 


3.26 


LaVOa 


m 


1 


0.81 


3.11 






3 


0.80 


3.12 



YBO3 compounds are characterized by smaller N° and, therefore, by somewhat larger U 
in comparison with the less distorted Lai^Os compounds. The same tendency holds for 
different transition-metal sites in the monoclinic structure: the sites with more distorted 
environment have large U and vice versa. 
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